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Abstract. The literature on vibrational thermodynamics of materials is reviewed. The emphasis 
is on metals and alloys, especially on the progress over the last decade in understanding differences 
in the vibrational entropy of different alloy phases and phase transformations. Some results on 
carbides, nitrides, oxides, hydrides and lithium-storage materials are also covered. 

Principles of harmonic phonons in alloys are organized into thermodynamic models for unmixing 
and ordering transformations on an Ising lattice, and extended for non-harmonic potentials. Owing 
to the high accuracy required for the phonon frequencies, quantitative predictions of vibrational 
entropy with analytical models prove elusive. Accurate tools for such calculations or measurements 
were challenging for many years, but are more accessible today. Ab-initio methods for calculating 
phonons in solids are summarized. The experimental techniques of calorimetry, inelastic neutron 
scattering, and inelastic x-ray scattering are explained with enough detail to show the issues of 
using these methods for investigations of vibrational thermodynamics. The explanations extend to 
methods of data analysis that affect the accuracy of thermodynamic information. 

It is sometimes possible to identify the structural and chemical origins of the differences in 
vibrational entropy of materials, and the number of these assessments is growing. There has been 
considerable progress in our understanding of the vibrational entropy of mixing in solid solutions, 
compound formation from pure elements, chemical unmixing of alloys, order-disorder transforma- 
tions, and martensitic transformations. Systematic trends are available for some of these phase 
transformations, although more examples are needed, and many results are less reliable at high 
temperatures. Nanostructures in materials can alter sufficiently the vibrational dynamics to affect 
thermodynamic stability. Internal stresses in polycrystals of anisotropic materials also contribute 
to the heat capacity. Lanthanides and actinides show a complex interplay of vibrational, electronic, 
and magnetic entropy, even at low temperatures. 

A "quasiharmonic model" is often used to extend the systematics of harmonic phonons to 
high temperatures by accounting for the effects of thermal expansion against a bulk modulus. 
Non-harmonic effects beyond the quasiharmonic approximation originate from the interactions of 
thermally-excited phonons with other phonons, or with the interactions of phonons with electronic 
excitations. In the classical high temperature limit, the adiabatic electron-phonon coupling can 
have a surprisingly large effect in metals when temperature causes significant changes in the elec- 
tron density near the Fermi level. There are useful similarities in how temperature, pressure, and 
composition alter the conduction electron screening and the interatomic force constants. Phonon- 
phonon "anharmonic" interactions arise from those non-harmonic parts of the interatomic potential 
that cannot be accounted for by the quasiharmonic model. Anharmonic shifts in phonon frequency 
with temperature can be substantial, but trends are not well understood. Anharmonic phonon 
damping does show systematic trends, however, at least for fee metals. 

Trends of vibrational entropy are often justified with atomic properties such as atomic size, 
electronegativity, electron-to-atom ratio, and mass. Since vibrational entropy originates at the 
level of electrons in solids, such rules of thumb prove no better than similar rules devised for 
trends in bonding and structure, and tend to be worse. Fortunately, the required tools for accurate 
experimental investigations of vibrational entropy have improved dramatically over the past few 
years, and the required ab-initio methods have become more accessible. Steady progress is expected 
for understanding the phenomena reviewed here, as investigations are performed with the new tools 
of experiment and theory, sometimes in integrated ways. 
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I. Principles and Methods 



1. Overview 

Without entropy to complement energy, thermodynamics would have the impact 
of one hand clapping. Boltzmann's graveside monument displays an equation for 
entropy of breathtaking generality 

S = k B lnn (1) 

(modernized slightly, with &b a s the Boltzmann constant). The nub of the problem is 
the number 0, which counts the ways of finding the internal coordinates of a system 
for thermodynamically-equivalent macroscopic states. Physical questions are, "What 
do we count, and how do we count them?" 

It is often possible to separate the internal coordinates of a material into con- 
figurational ones and dynamical ones. As an example, when the number Q for 
configurational coordinates counts the spatial arrangements of atoms over the sites of 
a crystal, this does not depend explicitly on temperature. On the other hand, when 
the number Q for dynamical coordinates counts the intervals of volume explored as 
atoms vibrate, this $7 increases with temperature, as does the vibrational entropy. 

Configurational entropies of atoms or spins undergo changes during chemical 
ordering or magnetic phase transformations, respectively. 1 Configurational entropy 
was largely understood by Gibbs, who presented some of the combinatoric calcu- 
lations of entropy that are used today [1]. The calculation of is more difficult 
when there are local correlations in the system, but cluster approximation methods 
have proved powerful and accurate [2, 3, 4]. An important point is that configura- 
tional coordinates have different equilibrium values at different temperatures, but 
temperature does not alter how the configurations are counted to obtain fi. 

Dynamical entropy grows larger as internal degrees of freedom of a solid, such as 
normal modes of vibration, become more strongly excited by thermal energy. 2 With 
increasing temperature more phonons are created, and the vibrational excursions 
of atomic nuclei are larger. Fundamentally, the entropy from dynamical sources 
increases with temperature because with stronger excitations of dynamical degrees 
of freedom, the system explores a larger volume in the phase space of position and 
momentum. This volume, normalized if necessary, is the Q for Eq. 1. 

Phonon thermodynamics, especially vibrational entropy, is the main focus of this 
review. It turns out that it is important to have rather accurate information on 
phonon energies or heat capacities to assess the importance of vibrational entropy. 3 

1 Electronic entropy can also have a configurational component in mixed- valent systems. Nuclear 
spins undergo ordering transitions at low temperatures, too, although at most temperatures of 
interest in materials physics the nuclear spins are fully disordered and their entropy does not 
change with temperature. 

2 Temperature also drives electronic excitations to unoccupied states, and when many states are 
available the electronic entropy is large. Spin excitations are another source of entropy, but care 
must be taken when counting them if the configurations of spin disorder are already counted. 

3 For example, if we apply Eq. 16 to a case where the phonon density of states curves of the 
a and (5 phases have the same shape, but differ in energy scaling by 10%, we obtain a change in 
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When comparing the phonon entropies of different alloy phases, a 1% accuracy in the 
difference in logarithmic-averaged phonon energy is desireable. This is challenging, 
so although there were many suggestions over many years that vibrational entropy 
could have a significant effect on the thermodynamics of order-disorder transforma- 
tions, until recently accurate computations of phonons were impractical, and the 
experimental results were ambiguous. Today, measurements and calculations have 
proved sufficiently reliable for useful assessments of the thermodynamic effects of 
vibrational entropy for different materials. This progress and the trends that have 
been found to date are the subjects of the present review. 

The independence of phonons is an important feature of the "harmonic model," 
for which there are many useful and powerful mathematical results. The harmonic 
model explains many physical phenomena, especially at low temperatures when the 
phonons are approximately independent of each other, and independent of other 
sources of dynamical entropy. It is often extended to higher temperatures as the 
"quasiharmonic model," which assumes independent phonons with altered frequen- 
cies. At high temperatures, however, the different sources of dynamical entropy can 
interact in important ways. The existence of phonons alters the energies required 
to create additional phonons (phonon-phonon interaction, PPI), and the electron 
excitations are altered by vibrational excitations (electron-phonon interaction, EPI). 
An understanding of vibrational entropy at high temperatures is emerging. 

The first sections of this review develop the principles of phonon thermodynamics 
and vibrational entropy. Methods for measurement and calculation are covered next, 
with emphasis on experimental methods. The third part is on the vibrational en- 
tropies of materials with different atomic configurations, and phase transformations, 
emphasizing low temperature behavior interpreted with a harmonic model. The 
fourth part is on non-harmonic behavior at higher temperatures, including the effects 
of electron-phonon and phonon-phonon interactions. 

2. Harmonic Lattice Dynamics 

The mass of an atom is in its point-like nucleus, and the forces between nuclei are 
transmitted through lightweight electrons, which are typically assumed to respond 
instantly to nuclear displacements (Born-Oppenheimer approximation). The Born- 
von Karman model assumes nuclear masses connected by massless electron springs, 
giving a set of ordinary second-order differential equations of motion. The model is 
specific to a crystal with translational periodicity, and its normal modes of vibration 
are wavelike and with dispersion relations between frequency and wavevector. The 
Born-von Karman model is expected to be reliable at low temperatures where atom 
displacements are small. 

Normal coordinates offer a sweeping simplification of the vibrational dynamics. 
In a harmonic solid, each normal mode does not transfer energy to the other modes. 
Without damping, a normal mode persists forever. This remains true in quantum 
mechanics, although integral numbers of excitations, phonons, exist in each normal 

phonon entropy of AS p h = 3&b ln(l.l) ~ 0.3fcB/atom. This change in vibrational entropy is almost 
half of the maximum possible change in configurational entropy of an order-disorder transformation 
of a binary alloy (which is fee In 2/ atom) . 
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mode. For harmonic solids, the energy cost for each additional phonon is £j = 
hwi {uji is the angular frequency of the oscillator), and is unchanged even as more 
phonons are added to each normal mode. Phonon modes can therefore be considered 
as independent harmonic oscillators, so we start with the thermodynamics of the 
quantum harmonic oscillator. 

2.1. Partition Function 

The partition function for a single harmonic oscillator of frequency Wj = Si/h is 

oo 

Zi = Ve^^ 1 /*' , (2) 



n 

e -/tei/2 



1 



-Pet 



(3) 



where Eq. 3 was obtained by identifying Eq. 2 as a geometric series times the 
constant factor exp(— f3ei/2), where (3 = (ksT) -1 . The partition function for a 
harmonic solid with N atoms and 3iV independent oscillators is the product of 
these individual oscillator partition functions, 



Zn 



3N -p £i /2 



-Psi 



(4) 



from which we can calculate the phonon free energy by the prescription F = 
-k B T In Z, 



-. 3N 3N 
^vib = ^}Z^ +fe B T^ln(l-e-^) , 



(5) 



and the phonon entropy (vibrational entropy) by differentiating with respect to T: 

3N 



S v[h = k Bz 2 -ln(l-e"^) + 



(6) 



Using the Planck distribution for phonon occupancy versus temperature, nfe, T), 

1 



n{e l ,T) = 



-In U 



-Pet 



e +PSi _ I ' 

l) = ln(l+n) , 
1 + n 



In 



n 



(7) 
(8) 
(9) 



it is straightforward to show that Eq. 6 reduces to the following form for the 
vibrational entropy per oscillator, S osc 



+/cb (1 + n) ln(l + n) — nlnn 



(10) 



Equation 10 can also be obtained from Eq. 1 with the combinatoric relationship 
SI = (M + m)\ I [Ml ml] for m phonons in M oscillators (with n = m/M). 
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Phonons are bosons, so it is interesting to compare Eq. 10 to the analogous 
equation for fermions, where m electrons are distributed over N states of energy e. 
For fermions, the combinatorics give f2 = N\ / [(TV — m)\ ml]. Fermion entropy is the 
same as the entropy of mixing of atoms on a crystal lattice, where each lattice site 
can hold only one atom. The entropy of mixing per atom is 



= —k B (1 — c) ln(l — c) + cine . 



(11) 



Equations 10 and 11 have two differences - all signs are reversed, and n can exceed 
1, whereas < c < 1. 

It is often useful to work with a phonon density of states (DOS), g(e). For a 3- 
dimensional solid, 3N 'g(e)de phonon modes are in an energy interval de. For a DOS 
acquired as digital data in m intervals of width Ae (so Ej = jAe), the partition 
function can be computed numerically 




3Ng( £j )Ae 



(12) 



which becomes at high temperatures 



m /, r r \3Ng(e j )Ae 



(13) 



Again, using Eqs. 7, 8, 9, we can obtain from Eq. 12 a useful expression for the 
phonon entropy of a harmonic material at any temperature (cf. Eq. 10) 



S vih (T) 



3£; B J #( £ ) 
o 



(e) + 1) In (n(e) + 1) - n(e) In (n(e)) 



de 



(14) 



where g(e) is normalized to 1 and n(e) is the Planck distribution (Eq. 7) for the 
temperature of interest. For thermodynamic calculations, series expansions of Eqs. 
5 and 6, valid at high temperatures, can be useful. For S p h(T), Fig. 1 compares the 
exact expression of Eq. 10 to approximations for up to four terms of a standard high- 
temperature expansion (terms are presented in the figure). At high temperatures, 
the following expression works well for the entropy of one oscillator mode with energy 
e = fko 



S vih ,UT) = k B (n(e) + l) In (n(e) + l) - n(e) In (n(e)) ~ k B 



,(15) 



although the constant factor of k B will cancel in differences of phonon entropy. 
A handy expression for the high-temperature limit of the difference in vibrational 
entropy between two harmonic phases, a and (3, can be obtained readily from Eq. 
14 



AS?- a = 3k B J (g a (e) - /(e)) ln(e) de . 



(16) 
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Figure 1. Temperature-dependence of entropy calculated for one mode of energy fko, with 
E = fc^T)- 1 . 

The important point about Eq. 12 is that the only material parameter relevant to 
the thermodynamic partition function in the harmonic approximation is the phonon 
DOS, g(e). Measuring or calculating the phonon DOS is central to much of our 
understanding of the vibrational entropy of materials. It turns out, though, that 
for most problems g(e) must be known with high accuracy, and this has been a 
challenge until relatively recently. 

2.2. Hamiltonian for Lattice Dynamics 

Since most of the atom mass is in the nucleus, for lattice dynamics we consider the 
Hamiltonian for nuclear motions in a crystal 

= + (17) 

l,K LTXlK 

and the characteristics of the atom vibrations require attention to the potential 
energy, <I>. Following the development of Maradudin, et al. [5] for a crystal lattice 
with a basis, the basis vectors are {r K ,n = 1, 2...T2.}. Each atom k, in the unit 
cell I, is allowed to vibrate about equilibrium with the displacement ui K {t). The 
instantaneous position Ri K (t) of atom Ik at time t is 

RUt) = n + f K + u lK {t) . (18) 

We use Cartesian components u a i K for the displacement vector, where cx — {x, y, z\. 
The total potential energy of the crystal, is a function of the instantaneous 
positions of all the atoms in the crystal, but it is expanded in a Taylor series of the 
atomic displacements about ui K (i) = (i.e., the equilibrium positions of the atoms 
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Rik = n + r K ) 

IkI'k' U a l K U a 'l' k' + ... , (19) 

aln a'l'n' 

where the coefficients of the Taylor series are the derivatives of the potential with 
respect to the displacements: 

<93> \ 

<$> a in = , (20) 

^aa' IkI'k' = q a ' ( 21 ) 

du a i K du a , VKl J q 

where the subscript zero means that derivatives are evaluated in the equilibrium 
configuration (all displacements equal to zero) and $o is the static potential en- 
ergy of the crystal. Because the force on any atom must vanish in the equilibrium 
configuration, we have [5] 

$«i« = V a,l,n. (22) 

In the harmonic approximation of lattice dynamics we keep only the remaining 
terms of the series written explicitly in Eq. 19 - we neglect terms of order three and 
higher in the displacements 

H n = ^ + ^0 + ^ E X! ®aa'lKl>K> U a l K U a 'l> K > , (23) 

Ik k alK a'l'K' 

where m K is the mass of the atom at basis index k of the unit cell. We rewrite the 
Hamiltonian in matrix form 

ff » = ^ + $ » + ^EE«t £wv > ( 24 ) 

Ik Ik V k? 

where a 3x3 force-constant sub-matrix is defined for each atom pair (Ik;1'k') 

QiKl'K' = [$* a >iKi'K>] ■ (25) 

If (1,k) / Eq. 21 applies. If (1,k) = (1',k'), ^aa'iniK is a "self-force 

constant," derived from the requirement of no overall translation of the crystal 

S/KiK = " E ^iKl'K' ■ ( 26 ) 

(J',«/)#(1.«0 

Because equal and opposite forces act between each atom of a pair, ^ K i> K i must be 
a real symmetric matrix. 
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2.3. Equations of Motion 

In the harmonic approximation, the equations of motion for all nuclei are 

m K ui K {t) = -^2 ®M'k' uv K '{t) V I, k . (27) 

1',k' 

In Eq. 27 there are 3 x TZ x N cel1 equations of motion to solve for a finite crystal 
containing N ccl1 unit cells. We seek solutions having the form of plane waves 4 of 
wavevector k, angular frequency u>^, and "polarization" e K j(k) (where j is a "branch 
index" discussed below) 



' Jl> e« J -(jfe)e i( *- fi - h 'fe t) . (28) 



2n(e?.,T) + l _ . lt ^ 



where we take the real part to obtain physical displacements, or we use the conven- 
tion of adding the complex conjugate of the exponential, and halving the prefactor. 
The phase factor, e lfc ' n , provides all the long-range spatial modulation of u lK ^.{t). 
The dependence on k, a short-range basis vector index, is taken out of the phase 
factor and placed in the complex constant e K j(k). It is convenient for the e K j(k) of 
Eq. 28 to have modulus unity, as does the exponential. The prefactors are consistent 
with the quantization of energy for one phonon, ftuj = m u> 2 (u 2 } = \ m u> 2 ti max , or for 
the thermal population of phonons in mode e^., having energy E = e^. \n(e^., T)+|]. 

(For root-mean-squared displacements, delete the y/2 in Eqs. 28 and 29.) 

We impose periodic boundary conditions on the finite crystal. These require the 
set of possible wavevectors {k} to have N ceYl values, a large number that gives a 
very fine mesh of fc-points in reciprocal space. With TZ atoms in the basis of the unit 
cell, the crystal has 3 x 1Z x N c ° n vibrational modes, in agreement with its total 
number of mechanical degrees of freedom. Each wavevector k is associated a-priori 
with 3 x 1Z types of vibrational modes, identified by a branch index, j. Each of the 
31Z different modes corresponds to a different polarization vector e K j{k) and angular 
frequency u^. (1 < j < 31Z), although degeneracies can be induced by symmetry. 

2.4. The Eigenvalue Problem for the Polarization Vector 

The polarization vector, e K j(k), is a characteristic of each vibrational mode k, j. The 
vector e K j(k) contains information on the excursion of each atom k in the unit cell 
for the phonon mode mode k,j. Specifically, it gives 1) the displacement direction 
of the atom, and 2) its phase lag in time with respect to the other atoms. The 
vectors e K j(k) for all the atoms in the basis (1 < k < 1Z) and their associated 



4 In general, the motion of an atom uu(t) will be a sum over many phonons u lt _^.{t) (cf., Eq. 
143). 
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angular frequencies w^. can be calculated by diagonalizing the "dynamical matrix" 

D(k). The dynamical matrix is obtained by substituting Eq. 28 into 27. It has the 
dimensions (3N x 3N) and is constructed from (3 x 3) submatrices D KK >(k) 



D(k) 



fD n (k) ... D 1N (k)\ 



\D N1 (k) ••• D NN (k)J 



(30) 



Each sub-matrix D_ KK i{k) is the Fourier transform of the force-constant matrix 
§Li K i> K i, considered as a function of (I' — I): 



D KK ,(k) 



1 



^/m K m K > ^7 



OkI'k' e 



ifc-(r ; ;-f()) 



(31) 



where we set I = because the summation is over all I' and the origin is arbitrary. 
By similarly collecting the polarization vectors into a vector of size 3 x 1Z, we rewrite 
the system of differential equations (Eq. 27) with the plane wave solutions (Eq. 28) 
as an eigenvalue problem: 



where 



D{k)e j {k)=u\.e j {k) , 
( e xlj (k)\ 

e ylj{k) 

e z ij(k) 
e X 2j(k) 

\e z Nj(k)J 



(32) 



ej{k) 



(33) 



It can be shown that the (31Z x 31Z) dynamical matrix D(k) is hermitian (for any 
value of k). It is fully diagonalizable, and the oo'i. are real. The 31Z eigenvectors and 

eigenvalues of the dynamical matrix evaluated at a particular wavevector k then 
correspond to the 31Z eigenmodes of vibration of the crystal for that wavevector. 



2.5. Calculation of the Phonon Density of States 

To calculate the phonon density of states (DOS) of the crystal, g(e), the dynamical 
matrix is diagonalized at a large number of points in reciprocal space (typically 
covering the first Brillouin zone). The diagonalization oiD{k) at each k point returns 
3TZ eigenvalues of angular frequency u^- (1 < j < 31Z), which are then binned into 
the DOS histogram. 

A phonon partial DOS, <7d(e), is a similar quantity, but it gives the spectral 
distribution of motion by one atom, the species d in the unit cell. Unlike the total 
DOS, g(e), the eigenvalues of the crystal are not weighted equally in the partial 
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DOS 

9d(e) = J2J2 6 d K \e aKj (k)\ 2 g(e) , (34) 
k aK i 

where the Kroneker delta is zero unless the atom species d is at the site k (when 
8dd = !)• Because the eigenvalues of the dynamical matrix are normalized for each 
k as 

J2K KJ (k)\ 2 = 1 , (35) 

the total DOS is the sum of the partial DOSs of all atoms in the unit cell, 

g(e)=Y.9d(e) ■ (36) 

d 

To calculate a partial DOS, for each diagonalization of the dynamical matrix at a 
specific k point, the partial DOS histogram for atom d is incremented at frequency 
by the amount J2 a e adj (^) e adj (k)- The partial DOS for atom d is large at energies 
where there are many modes with large displacements of atom d. The histogram for 
the phonon partial DOSs are computed simultaneously with the total DOS. 



3. Predictions with the Harmonic Model 

The harmonic model provides much of the truth about vibrational entropy at low 
temperatures, and its analytical formulation allows a number of important pre- 
dictions about vibrational entropies of alloying or phase transformations. Some 
important results are presented here without full proof. Two classic and sophis- 
ticated treatments of the harmonic solid are Maradudin, et al. [5], and Born and 
Huang [6]. 

3.1. Mass, Local Modes, and Resonance Modes 

Otto Stern reported a number of theoretical observations relevant to lattice dynamics 
in 1916 [7]. Importantly, he showed that the trace of the dynamical matrix is equal 
to the product of the normal mode frequencies squared. From this it is possible to 
show that mass disorder does not contribute to the vibrational entropy of order- 
disorder or unmixing transformations, for example. If a new structure is formed by 
rearranging the same atoms, the effects of the masses on the vibrational entropy can 
be neglected, at least in the classical limit. All changes in vibrational entropy are 
attributed to changes in force constants (or changes in the interatomic potential), 
and not to rearranging the same masses [8, 9, 10, 11] 

The kinetic energy for each atom is p 2 /(2m), where the momentum, p, is a 
maximum when the atom is at its equilibrium position (ui K = 0) and all the thermal 
energy is kinetic. Rearranging the different atom masses into different structures 
has no effect on the kinetic energy of each atom as it moves across its equilibrium 
position. There is no change in the volume of momentum space explored during 
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vibration. The potential energy of each oscillating atom does not depend on mass, 
only position. Changing the mass of an atom will not change the range of its motion 
in the potential well (Max[u/ K (i)]), or the volume of phase space in the position 
coordinates. Since the volume of momentum-position phase space is unchanged by 
rearranging the same masses, by Eq. 1 there is no change in vibrational entropy 
in the classical limit. Only changes in interatomic force constants are important 
when comparing classical vibrational entropies of different phases in the harmonic 
model. Anharmonicity does not change this independence of mass, but it alters other 
results. 5 

The Einstein model, which assumes each atom is an independent harmonic oscil- 
lator, cannot account for wavevector-energy relationships of phonons. Nevertheless, 
the Einstein model can be useful when little energy is transferred between the atoms 
in the solid, and there is little change of phonon energy with wavevector. Consider 
hydrogen atoms in interstitial sites between metal atoms. Owing to its low mass, a 
hydrogen atom vibrates as if in a cage. Although the neighboring metal atoms must 
move a bit in response to the moving hydrogen atom, these metal atom displacements 
are small. If the interatomic forces between the metal and hydrogen atoms arc 
moderate, the system will show a "local mode" at a high frequency, typically above 
the modes of the metal crystal. The local mode is dominated by motions of the 
hydrogen atoms, although the metal atoms must have some contribution in their 
own partial density of states. 

The other extreme of a heavy impurity atom in a host of lighter atoms produces 
a "resonance mode," where the partial phonon DOS of the heavy impurity is dom- 
inated by a peaked spectrum at low energy. The neighboring matrix atoms move 
in response to the local motions of the heavy impurity atom, but the dynamics 
is not simple. Because the resonance mode exists at the same frequency as other 
modes in the matrix, there is substantial coupling between the resonance mode and 
the continuum, so the resonance mode is not an eigenmode of a harmonic solid. A 
resonance mode appears as being highly damped, with a short lifetime and a broad 
width in energy, sometimes of the form [5] 



The resonance mode is not intuitive, and its discovery was an achievement of 
theoretical work [12, 13, 14, 15]. 

The existence of a resonance mode depends on the ratios of the masses and 
interatomic force constants. If all interatomic force constants are equal, a resonance 
mode will exist in an fee lattice if the impurity mass is about 2.6 times heavier than 
the atoms in the matrix [5]. More likely, the interatomic force constants will differ 
between the host atoms and the impurity, compared to host-host force constants. 
Stiffer forces to the neighbors of the impurity effectively increase the mass of the 
defect center, favoring a resonance mode. 

° For harmonic oscillators, the thermal energy is equipartitioned between kinetic and potential 
energy. The equipartition of Icb T/2 into each degree of freedom of momentum and position origi- 
nates with these coordinates appearing quadratically in the Hamiltonian H, producing a Gaussian 
integration of e ~ H / k B T \ A cubic or quartic term in H causes an imbalance between kinetic and 
potential energy. 




(37) 
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3.2. Long Wavelength Limit and the Debye Model 

The relative motions of different atoms in the unit cell are easiest to understand in 
the long wavelength limit. Each inequivalent atom in a unit cell forms a sublattice 
connecting its equivalent atoms in all other unit cells, distinct from the sublattices 
of the other atoms. In the limit of long waves (as the wavevector k — > 0), these 
sublattices move rigidly, as seen on the scale of many unit cells. In optical modes, 
the sublattices move against each other, out-of-phase for example, but equivalent 
atoms in adjacent unit cells have the same displacements at the same time. In the 
acoustic modes, the motions of all sublattices are the same, as adjacent unit cells 
distort their shapes in unison. 

Acoustic modes at long wavelengths are elastic waves through solids. The anisotropy 
of a crystalline solid adds complexity to the problem, but consider the wave equation 
for the propagation of a longitudinal wave along the [100] direction: 



Relationships like Eqs. 38 and 39 hold true in all directions, although different combi- 
nations of elastic constants are required. Equation 39 requires a set of relationships 
between the elastic constants (such as Cu) and the interatomic force constants 
^ikl'k'i which determine the phonon dispersions through Eqs. 31 and 32. Physically, 
the elastic constants originate from the interatomic forces, so such a relationship is 
expected. The derivations of such relationships are not simple, however, owing to 
the anisotropy of the interatomic forces and differences between crystal structures. 

Sound velocities are the rigorous limit of dui/dk as k — > 0. At the atomic scale, 
a long wavelength longitudinal phonon is a density wave of compression followed 
by tension. In a local region (much smaller than A), all atoms move collectively 
towards and against their neighbors, irrespective of their mass or interatomic force 
constants. When A exceeds the range of the interatomic forces, the wave propagation 
becomes independent of atomic detail, and in a given direction all waves of the same 
polarization propagate with the same speed. Continuum behavior is achieved. 

A subtlety arises when anharmonicity is significant, such as at high temperatures. 
Anharmonicity can suppress phonon lifetimes to tens of vibrational periods or less 
(Sect. 19.5). In this case, the long acoustic waves measured by ultrasonics propa- 
gate through a medium influenced by averaged anharmonic effects. On the other 
hand, the acoustic phonons measured by neutron inelastic scattering propagate in 
their reduced lifetimes across a medium where equilibration of vibrational energy 
is not taking place between phonon modes. With anharmonicity, the energies and 
du/dk of the neutron measurements are often larger than expected from ultrasonic 
measurements. 

The Debye model of phonons in solids extrapolates this long wavelength behavior 
to all phonon frequencies, to, by assuming they are linear functions of wavevector, k = 




which has a wave speed of vuioo] = \/C\\j p {C\\ is an elastic constant and p is 
density). This wave speed i^ioo] must equal the group velocity obtained from the 
longitudinal branch of the phonon dispersion along [100] at small k: 



Lira dw/dfc 

fc [ioo]-*0 




(39) 
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2tt/X. This assumption never works for optical phonons, and is erroneous for acoustic 
waves when k approaches a reciprocal lattice vector (i.e., A is interatomic distances). 
In particular, the highest phonon frequency occurs when nearest-neighbor atoms 
move out-of-phase with respect to each other. Here there is no wave propagation, 
and doo/dk = 0. 

Comparing the Debye temperatures of phases before and after a phase transfor- 
mation allows an estimate of the change in vibrational entropy, using for example a 
high temperature expression like AS 01 ^ 13 = 3A?b ln(0Da/©D/3), where Q^a and ©d/3 
are Debye temperatures of phases a and (5. This must be done with care, however, 
because the Debye model has a history of misleading the materials physics commu- 
nity about the importance of vibrational entropy. The Debye model is unreliable for 
any calculations of vibrational entropy when the Debye temperature is obtained from 
experimental measurements of long wavelength behavior such as sound velocities, 
because only the low-energy tail of the phonon DOS is sampled. Differences in 
this tail usually do not correlate reliably with differences in the phonon spectra 
of different materials because most of the phonons have wavelengths of atomic 
distances. On the other hand, Debye temperatures obtained from heat capacity 
measurements are more representative of the phonon contributions to the entropy, 
and are more trustworthy for estimates of vibrational entropy. Beware that the 
Debye temperature is sometimes graphed as a function of temperature, which is 
a confession that the high frequency phonons cannot be obtained reliably from 
long- wavelength vibrations. 

Fernandez Guillermet and Grimvall tested modifications of the Debye tempera- 
ture for evaluating vibrational entropy and phonon spectra [16]. Specifically they 
defined a Debye temperature from the standard entropy of pure elements, as ob- 
tained for example from the JANAF Thermo chemical Tables, with a small correction 
for electronic entropy. A second Debye temperature was obtained from bulk moduli, 
but it differed substantially from that for the heat capacity. For Ru, for example, the 
two Debye temperatures were 384 and 555 K, and the effective interatomic force con- 
stants derived from the two Debye temperatures differed by a factor of two. 6 For M, 
4d, and 5d elements in the same column of the periodic table, there are fairly constant 
ratios of the two types of force constants. Nevertheless, the large disparities between 
Debye temperatures obtained from the entropy and from the elastic constants do not 
bode well for determining vibrational entropy from Debye temperatures. For some 
bcc and hep transition metals, the ratios of Debye temperatures deduced from the 
entropy at ambient temperature and at cryogenic temperature range from 0.8 to 1.0 
[17]. This indicates that the Debye temperature useful for determining differences in 
vibrational entropy is obtained from the vibrational entropy at the temperature of 
interest, i.e., the Debye temperature should be a recasting of the vibrational entropy 
itself. 

It is possible to account for the effects of thermal expansion on the Debye temper- 
ature with a Debye-Griineisen approximation, where Qd(T) = 0^ \ro/r(T)] 37 . With 
increasing T, 0d(T) becomes smaller as r(T) becomes larger. This approach was 
used by Moruzzi, et al. to calculate thermal properties of elements [18]. Sanchez, et 
al., used this approach to improve the effective cluster interactions for first-principles 

6 The force constant, k — Mlu 2 , and with tko = fc B 6D, k = M(k B Q-D/K} 2 . 
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phase diagram calculations [19]. Tseng and Stark used this Debye-Griineisen approx- 
imation to calculate the solubility of Ti in Al at both elevated temperatures, and 
under pressure [20]. 

3.3. Disordered Systems 

The pioneering analytical work on the disordered two-component chain by Dyson 
[21] inspired a number of studies on the dynamics of disordered systems. This field 
is excellent for computational methods, and the value of computer calculations of 
lattice dynamics was realized as soon as appropriate computing resources were 
available. Early numerical computations on disordered systems were reviewed in 
1972 by Dean [22] in a computational science article that is still interesting today. 
Much later, a pedagogical computational study of phonon propagation in disordered 
chains was performed by Allen and Kelner [23]. They showed that with increasing 
disorder, phonon propagation evolves from ballistic, to diffusive, to localization of 
the vibrational energy. The transitions between these three regimes are not abrupt. 
With modest disorder, some phonon wavepackets that propagate ballistically can be 
seen to scatter off the stronger inhomogeneities, and again propagate ballistically 
to their next scattering. These computational studies, although informative in their 
own right, were not directed towards the study of vibrational entropy. 

4. Bond Proportion Model 

If two phases differ in their energies by AE and entropies by AS, these two phases 
will have the same free energy (i.e., AF = AE — T C AS = 0) at the critical 
temperature T c , determined by the condition AE = T C AS. Now suppose that T c 
has been evaluated for the difference in configurational entropy alone. Calling it 
T c conf , we have AE = T c conf A5 conf . Now if, for the same AE, we add a difference 
in vibrational entropy AS^ib to the model, the critical temperature of the phase 
transformation is changed to T£° +Vlb . j n classical limit, this new critical tem- 
perature is related to the critical temperature with configurational entropy alone as: 
T conf+vib (AS^f + AS vih ) = T™ nl AS coni . There are handy methods for evaluating 
AE and A5 con f for configurations of atoms on a lattice. A challenge is to determine 
how ASVib depends on the configurations of atoms, so that a new critical tempera- 
ture can be obtained for a more complete model of a phase transformation of atoms 
on a lattice. 

The vibrational entropy of phase transformations depends on changes in bond 
stiffnesses (force constants), so a natural first approach to estimate the vibrational 
entropy of a phase transformation is to count the different types of bonds between 
neighboring atoms, and how they change after a phase transformation. With a little 
risk, the accounting can be confined to first-nearest-neighbors (INN), and with more 
risk we can associate the stronger bonds with the stiffer springs. This reasoning is 
enough to explain qualitatively the effects of vibrational entropy on the critical tem- 
perature, where a disordered phase at high temperatures transforms upon cooling 
to either an unmixed or ordered structure. The disordered phase has more bonds 
that are energetically unfavorable, and therefore weaker springs between atoms. A 
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disordered phase is expected to have the larger vibrational entropy, increasing its 
relative stability. The critical temperature for disordering is therefore expected to 
decrease when vibrational entropy is added to the configurational entropy, and this 
is the usual case for real alloys. Assigning a fixed stiffness to a particular type of 
chemical bond, and counting these bonds, is called the "bond proportion model." It 
has a long history [9, 24, 25, 26, 27, 28]. 

Unfortunately, for quantitative predictions of the vibrational entropies of alloy 
phases, the bond proportion model is rarely adequate. Local atomic configurations 
alter the bond stiffnesses of atom pairs, and Section 5 discusses such effects parame- 
terized by bond lengths and materials chemistry. Chemical effects are also discussed 
in Sections 18-20, including effects of temperature. 

4.1. Bond Proportion Model and Unmixing on the Ising Lattice 

Consider a binary alloy with atom species A and B arranged randomly on the 
sites of a crystal. The concentration of B-atoms is c, and the concentration of A- 
atoms is 1 — c. The crystal has N sites, and therefore a number cN of B-atoms and 
(1 — c)N of A-atoms. Each site is surrounded by z sites as first-nearest-neighbors 
(INN). For simplicity, we make the "point approximation," and assume that the 
chemical potential of an atom on a site depends on only the average occupancy of 
its neighboring sites. For a specified overall composition, we obtain thermodynamic 
functions with only knowledge of the species of atom on each site. 

To obtain the total energy of an alloy configuration, we count the number of INN 
pairs: A-A, B-B, (A-B = B-A), each with energy e AA , cbb> e A B- For a random solid 
solution, on average each A-atom is surrounded by zc B-atoms, and its other (1 — c)z 
neighbors are A-atoms. Attributing all of the energy of the average A-atom, e A , to 
its INN chemical environment 

z z 
e A = -(1 - c) e AA + -ce AB , (40) 

Z Z 

e B = ^ ~ c ) eAB + 2 CeBB ' ( 41 ) 

(and likewise for the average B-atom energy, The division by 2 corrects for 
the double-counting that occurs when performing the average by summing over all 
atoms in the crystal (because the bonds from the neighbor atoms back to our central 
atom are counted again) . The energies e A and e-g also include the vibrational energy 
of the A- and B-atoms, which is 3k^T in the classical limit. Since this is the same 
for both types of atoms, however, it does not affect the present analysis of phase 
transformations on an Ising lattice. 
We define the "interchange energy" 

V = (e A A + e B B -2e AB )/4. (42) 

The physical meaning of V is made clear when we consider the interchange of an A- 
B pair as shown in Fig. 2. For the initial configuration of the central A and B-atoms 
in Fig. 2a, we count 5 A-B pairs, 1 A-A pair, and 1 B-B pair. After interchange of 
only the central A and B-atoms, 7 we have 3 A-B pairs, 2 A-A pairs, and 2 B-B pairs 



Note that this interchange conserves the total number of atoms of each type. 
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Figure 2. (a) Initial configuration around central A and B-atoms. (b) New configuration after 
exchange of the central A-B pair. 



(Fig. 2b). The change in energy is precisely W = eAA + e BB — 2eAB (or an integer 
multiple of 4V for different configurations). 

The partition function for a single site, ^i s ite is 



'lsite 



e -/3e A + e -/3e B j (43) 



where (3 = l/{k^T). Assuming all site occupancies are independent, the total 
partition function of the alloy with N sites is 

Z N = Z^ e . (44) 

We can evaluate Zn by the binomial expansion of Eq. 44 plus Eq. 13 for the high 
temperature limit. For a fixed n = cN, we select only one term of the binomial 
expansion: 



3N g(tkjj) HAu) 



J N 



(N - n)\n\ 



Evaluating the last factor in Eq. 45, the partition function for the vibrational fre- 
quencies, requires an approximation. We expect a physical correspondence between 
nearest-neighbor pairs and interatomic force constants. Section 2.4 showed how the 
normal mode frequencies {uJj} are obtained from coupled equations of motion, and 
in principle a full diagonalization of the dynamical matrix may be required for each 
composition. This is difficult in practice, since disordered solid solutions cannot be 
described as periodic structures with small unit cells. On the other hand, we note 
that in Eq. 45 we need products of all the frequencies, not necessarly the individual 
frequencies. Mahanty and Sachdev [29] and Bakker [30, 31] addressed this problem, 
and showed for some crystals how the product of the equals the product of all 
interatomic force constants, {7^} between atom pairs, i — j, divided by the product 
of all atom masses, which is constant in a closed system. This result is consistent 
with the following approach of adapting the Einstein model for vibrational dynamics 
to our thermodynamic point approximation. 

In the Einstein model each atom vibrates independently, as if tethered to a lattice 
site. 8 Its motions in three dimensions are decoupled, so each atom mass has three 
independent springs. The correspondence works naturally for a linear chain, square 



The mass of the oscillator is therefore the mass of the atom, not the reduced mass of the atom 
plus its neighbor. In this approximation, the Einstein model preserves the independence of mass 
discussed in Section 3.1. 
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lattice, or simple cubic lattice, for which each atom is associated with 1, 2, or 3 
interatomic springs, with a frequency Wjj that depends on the INN pair of atoms, i 
and j. For a simple cubic lattice with z = 6, each atom has z/2 neighbors (bonds are 
not overcounted), and z/2 springs (normal modes are not overcounted). When the 
coordination number z is not the same as 2V, where V is the number of dimensions, 
we correct the mode counting with the ratio rj = 2D/z, which is 1/2 for an fee crystal, 
for example. The vibrational factors in Eq. 45 now include only the frequencies 
{uj aa , ^>bb, ^ab}- Their exponents for the numbers of modes are replaced by the 
numbers of bonds of types A-B, A-A, B-B, adjusted to conserve modes, i.e., rjN AB , 
r)N AA , r]N B B, giving 



Z N (c) 



AM 



-/3e A 



N-n 



( e -/3e B y 



(N - n)\n\ 

( k B T y NAB / k B r y NAA / k B T y N ^ 

\haj AB ) \hw AA ) \hw BB J 



(46) 



To obtain A^aa for a random solid solution, we consider the average number of 
A-atom neighbors of an A-atom to be (1 — c)z, so: 



N AA = (l-cf Z 4 



(47) 



where dividing by 2 avoids double-counting. Likewise, Abb = c 2 zN/2 and Aab = 
c(l — c)zN (where we add A-B and B-A bonds together). 

Obtaining the Helmholtz free energy from the expression F(c) = — k B TlnZ(c), 
and using the Stirling approximation ln(x!) ~ xlnx — x, we obtain from Eq. 46 



where: 

-Econf(c) 

'S'conf ( c ) 
5"vib(c) 

'S'vib(c) 



F(c) = E con{ (c) - T S'conf (c) + Svib(c) , 



e A A + 2c{e A B ~ e A A) + c 2 {AV) 
-k B \clnc+ (1 - c) ln(l - c)]n , 



zN 



k B C AA + 2c{C AB - C AA ) + c 2 (C AA + C BB - 2C AB ) 



ZT]N 



k B 



ln^ + 2cln^ + c 2 ln^i^ 



hv AA 



^AB 



WAA^BB 



VN. 



(48) 

(49) 
(50) 
(51) 
(52) 



Here i? con f(c), S'conf (c) and S v ;b(c) are the configurational energy, configurational 
entropy, and vibrational entropy of the random alloy. In Eq. 51 factors were defined 
such as C AB = ln(k B T/(hu> AB )). In the classical limit of the Planck distribution, 
n(e AB , T) of Eq. 7, this C AB is the logarithm of the number of phonons of frequency 
u> AB . Additionally, k B C AB is the vibrational entropy of an oscillator of frequency 
waBj Eq. 10 in the classical limit. Equations 49 and 51 show that E con f(c) and 
S v ib(c) have identical quadratic forms in c, reflecting the energies or interatomic force 
constants of nearest-neighbor bonds in the point (Bragg-Williams) approximation. 

For chemical unmixing we consider the case where V < (from Eq. 42 this 
means that the average of A-A and B-B pairs is more favorable, i.e., their energies 
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Figure 3. (a) 5 , CO nflg(c) for random solid solution (Eq. 50), -Econfig(c) for unmixing, i.e., V < (Eq. 
49), representative 5 , v ib(c) when lj% b < uaa^bb (Eq. 52), and 5 to t(c) = 5 CO nHg(c) + 5 v ib(c). (b) 
Free energy for zV/(kT) — 1 (critical condition at c = 0.5 without vibrational entropy), showing 
effects of different ratios of cj%b to ojaa^bb- 



are more negative, than A-B pairs). 9 The energy and entropy terms of Eqs. 49, 
50 and 52 are graphed in Fig. 3a. It was assumed that cjab = 0.77- v /u;aawbb, arid 
ujaa < ^bb, as may be plausible if eAA > eeB- Figure 3b shows curves of the 
total free energy versus composition for different ratios of vibrational frequencies. 
The logarithmic singularity of «S , con fig(c) always forces a downwards slope away from 
c = and c = 1. The energy term, positive in this case with V < 0, promotes 
unmixing at intermediate concentrations at low temperatures. 

The critical temperature for chemical unmixing, T c , is found by setting d 2 F/dc 2 = 
0. It is straightforward to do this with Eqs. 49 through 52, with the result 



k B T c 



A\V\z 



_ [ _L_ + 2 P ln(^M^ 



"AB 



(53) 



The T c is reduced by vibrational entropy when waa^bb > <^ab > anc ^ increased when 



^AA^BB < ^AB- 



When 



AB 



cjaa^bB) Eq. 53 provides T c = z\V\/kB for c = 0.5; the familiar 



result without vibrational entropy. The middle curve in Fig. 3b shows F(c) for this 
condition. The free energy without vibrational entropy, F con fi g (c), has an interesting 
curvature in c below T c = z\V\/kB, which can be helpful for understanding spinodal 
decomposition. 



9 The opposite sign of V > means that at low temperatures the alloy will tend to maximize 
the number of its A-B pairs by developing chemical order, as described in the next section. 
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When w AB / cjaa^bb, <Svib(c) is a quadratic function of c (Eq. 52), much like 

-^config (c) • 

When u\ B < ljaa^bB) the vibrational entropy term — T5 v ;b(c) effectively dimin- 
ishes -Econfig(c), reducing T c . 

When u\ B > ojaa^bB; T c is raised and the whole unmixing phase bound- 
ary moves to higher temperatures. If w AB exp(— 2/3) > cjaa^bb (approximately 
0.513 u\ B > waa^bb)) disordering cannot occur at any temperature at c = 0.5. 
For yet larger w AB with respect to ljaa^bb, the entropy of mixing quickly becomes 
negative for intermediate compositions, suppressing disordering. A negative entropy 
of mixing has been observed in alloys of V-Pt at surprisingly low concentrations of 
Pt (see Sect. 20.3.2). 

4.2. Bond Proportion Model and Ordering on the Ising Lattice 

For chemical ordering, we consider the case V > of Eq. 42, where the energetics 
tends to maximize the number of A-B pairs at low temperatures. At high tempera- 
tures, we expect this alloy to become a random solid solution if the configurational 
entropy is dominant. For illustration of an ordering transformation we choose an 
equiatomic A-B alloy on the square lattice shown in Fig. 4. Figure 4 identifies the 
a- and /3-sublattices, each having A/2 sites for a crystal with N sites. Finite tem- 
peratures favor some disorder in the structure - putting a few atoms on the wrong 
sublattice gives a big increase in the configurational entropy, 11 driving some mixing 
of A-atoms onto the /3-sublattice, in exchange for B-atoms on the a-sublattice. 

Figure 4 motivates the definition of four types of sublattice concentration vari- 
ables: {A A , A A , TVg, Ag}, where for example A A is the number of A-atoms on the 
a sublattice. For an equiatomic alloy, by conservation of atoms and conservation 
of sublattice sites, N£ = N B = R, where R denotes atoms that are on the "right" 
sublattice, and A A = N B = W, denoting "wrong" atoms. The sublattice concen- 
trations are used to define a "long-range order" (LRO) parameter, L, for the alloy 
[32, 33, 34] 

This L can range from — 1 < L < 1. The case of L = 1 corresponds to perfect order 
where R = N/2 and W = 0. The disordered alloy has as many right as wrong atoms 
on each sublattice, i.e., R = W, so L = for a disordered alloy. 

Analogous to the assumption of a random solid solution, we assume that the 
atoms are placed randomly on a sublattice, e.g., the A-atoms have no preference 
for any site on the a-sublattice, irrespective of the presence of A-atoms or B-atoms 
on neighboring sites. This is the "point" or "Bragg-Williams" approximation. The 
strength of this approximation is that it makes it easy to count the bonds between 
the different species of atoms. 



10 A simple approximation is that the strengths of the bonds and spring constants are propor- 
tional, giving •S'vibfc) the same shape as -E C onfig(c). This is conceptually convenient, and drove some 
early studies, but it has weak predictive power. 

11 The entropy has a logarithmic singularity in sublattice concentrations analogous to the 
concentration dependence of Eq. 50 near c = 0. 
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Figure 4- Ordered square lattice with two species {A, B}. Dashed and solid lines represent two 
interpenetrating square sublattices {a, f3}. A few misplaced A and B atoms are shown. 



To determine the numbers of pairs A^aa, -Abb, Nab, an important point evident 
from Fig. 4 is that the first neighbors of each atom on the a-sublattice are on the 
/3-sublattice, and vice-versa. When evaluating -/Vaa we start with all R A-atoms on 
the a sublattice, and multiply by the probability of finding a wrong A-atom in each 
of the z neighboring sites of the (5 sublattice (which has N/ 2 sites) . 



N A A = N BB = RW 



N/2 ' 



A^AB = (R 2 + W 2 



Using Eqs. 54 through 56 



A^aa = N BB = (1 - L 2 ) 
N AB = (l + L 2 ) 



N/2 



zN 



zN 



(55) 
(56) 

(57) 
(58) 



In the high temperature limit for phonons, we write the partition function for the 
alloy in terms of the random occupancies on two sublattices (giving two binomial 
probabilities) and three factors that depend on the numbers of bonds: 



f (N/2)\ \ 



p e AA fMLJ 1 



N A A 



-P CAB 



huJAA 
huJAB 



p e BB (MLY 

huJBB/ 



N B b 



(59) 



As for Eq. 46, we seek a correspondence between the number of modes and the 
number of pairs. For the square lattice Eq. 59 works naturally with rf = 1, but 
for other lattices it is necessary to take care to avoid incorrect counting of the 3N 
normal modes by setting i] = 2D/z, where T> is the dimensionality. 

The partition function depends on the LRO parameter, L, and we replace the 
exponents in Eq. 59 with Eqs. 57-58, and use Eq. 54 with the conservation of sites 
(R + W = N/2) to remove all instances of R, W, N A a, N B b, N A b in Eq. 59, 
obtaining Z(L). The free energy is obtained as F(L) = — ksT In Z(L). Subtracting 
the free energy of the disordered solid solution, we obtain the free energy of ordering, 



20 



B. Pultz 



AF(L) = F(L) - F(L = 0) 

v N 

AF(L) = -L 2 —4V 



-T jk B [-(1 + L) ln(l + L) - (1 - L) ln(l - L)] 



2 VN r, i„ ( ^AA^BB^ 



-TL 2 ^- k B In A T , (60) 



AB 



I 3 



where the three lines of Eq. 60 are identified as the changes upon ordering of 
-^config) —TS con fig, and — TS^, respectively. We set dAF(L)/dL = and obtain 
a transcendental equation for the equilibrium value of L(T) 

2LzV = k B T\n(l±4) + LVk B T\J ^ B ). (61) 
\l-LJ \uaau BB J 

This differs from the usual Bragg- Williams solution by the extra term from vibra- 
tional entropy. This term is linear in L, so it can be considered as a temperature- 
dependent rescaling of the interatomic exchange energy, V. 

There is no analytical solution to Eq. 61 for all T, but it is possible to obtain 
an analytic solution for the critical temperature, T c . The trick is to assume L is 
infinitesimally small, so Eq. 61 can be linearized. 12 The result is 

k B T c = Z -^- 2 . (62) 

In the absence of vibrational entropy, i.e., u\ B = ljaa^bb, Eq. 62 reduces to 
T c = zV/k B , a well-known result of Bragg- Williams theory. The configurational 
entropy is always important in a chemical order-disorder transformation. When 
both vibrational and configurational entropies are important, if uj\ b > waa^bb 
the Bragg- Williams value of T c is reduced, perhaps significantly. In the case where 
oj\ b = I.Iwaa^bB) f° r a simple cubic lattice T c is reduced to 87% of its value 
without vibrational entropy. This suppression of T c seems to be the usual effect 
of vibrational entropy, but the opposite is possible when lo\ b < waa^bb- Al- 
though difficult to expect from physical intuition, in the case where w^g/cJAA^BB < 
exp(— 2/3), vibrational entropy suppresses disordering, and T c = oo. 



4.3. Monte Carlo Results 



Some of the earliest computational studies on the vibrational entropy of ordering 
were performed by H. Bakker [30, 31]. Building on work by Mahanty and Sachdev 
[29] , Bakker used the result for the difference in vibrational entropy of phases a and 

12 This assumes that the ordering transformation remains second order, which is true for Eq. 
61. The curve of L(T) from Bragg- Williams theory, which decreases continuously from L — 1 to 
L — at T c , is rescaled monotonically with SVib - the effective V can be considered as altered by 
a temperature-dependent term. 
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ASUr = 1/2 In ^ 




(63) 



where the determinant of the mass matrix M_ is the product of the atomic masses, 
and undergoes no change with ordering. The force constant matrix is 5. In the case 
of a simple cubic A-B alloy with first-neighbor interatomic forces 7AA, 7bb, TAB, 
the resulting critical temperature is 



in agreement with Eq. 62 because 7 ~ to 2 M. 

Tuijn and Bakker [35, 36] computed the effects of vibrational entropy on the 
heat capacity of a three-dimensional crystal by counting the number of A-B pairs 
in a Monte Carlo simulation. Tuijn's results are shown in Fig. 5 for various values 
of A = 7ab/[7AA7bb] -1 / 2 - As expected, with A > 1, and hence stiffer vibrations 
in the ordered phase, the critical temperature is reduced. The critical temperature 
is reduced by approximately a factor of 1.8 when A = 2. The transformation still 
appears to be continuous, i.e., second order, although the singularity is rounded by 
the finite size of the lattice used in the simulations. Figure 5 shows that the more 
accurate Monte Carlo simulations predict a lower T c than Eq. 64, a well-known 
result without vibrational entropy. The relative sensitivity to A is also different for 
the Bragg- Williams and Monte Carlo results, however. The best agreement between 
a rescaled Eq. 64 and Fig. 5 is found for z = 9, rather than 6. Tuijn and Bakker also 
explored the effects of second nearest-neighbor interatomic forces and non-central 
forces. The heat capacities were similar to those of the simple cubic lattice with 
first-neighbor interactions, except when A>2. 



5.1. Phonon Frequencies and Bond Lengths 

The bond proportion model just described, which counts the different types of atom 
pairs, has a deficiency that is easy to understand. The bond lengths between pairs of 
atoms are different in different structures, and the stiffnesses of a bond is sensitive 
to the length of the bond. 13 Counting the proportions of atom pairs is important 
when comparing the vibrational entropies of two different structures, but it is also 
important to know how each bond stiffness changes with bond length. 

Consider a disordered solid solution with a distribution of interatomic distances, 
shorter and longer than optimal, owing to the different atomic sizes. The calculation 
of these distances is not simple because the relaxations in position for each atom 
depend in turn on the relaxations of its neighbors. Nevertheless, we can estimate 



zV 



(64) 




5. Bond-Stiffness-versus-Bond-Length Model 



1 This same characteristic of the interatomic potential also gives rise to quasiharmonic theory, 
where interatomic forces weaken with thermal expansion. 
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-T/v (l0 4 K/eV) 



Figure 5. Heat capacity at constant volume for harmonic solids with different ratios of A, as a 
function of —T/v = T/2V in the present notation. B.W. denotes Bragg-Williams approximation. 
After [37]. 



the size of this effect on phonon frequencies. The Griineisen parameter 7 is defined 
as: 

^~Zdv> (65) 

where to is the phonon frequency and V is volume. Recognizing that V oc r 3 and 
uj oc \fkl, where k s is an interatomic force constant 

r dk s , . 

«7 - ~~ ■ (66) 

There is no universal value for a Griineisen parameter, but typically 7 ~ 2. In 
this case, Eq. 66 predicts that if the interatomic distance changes by +1%, the 
interatomic force changes by -12%. This is a major effect. A change in interatomic 
distance of +10% could cause the interatomic force to be lost entirely. 

Grimvall and co-workers [38] correlated the interatomic forces of 8 alkali halides 
with the NaCl structure to their volume per atom in the crystal structure. (The 
interatomic forces were deduced from the bulk moduli.) The interatomic forces 
decreased approximately inversely with atomic volume per atom as k s = 1300/F 
where the force (spring) constant k s has units of N/m and V has units of 10 -30 m 3 . 
Fernandez Guillermet and Grimvall sought similar trends in vanadium and titanium 
compounds with oxygen, nitrogen, and carbon [39]. The trend in this case was a 
rapid increase in the k s with V, and there was a significantly different slope for the 
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titanium and vanadium compounds. Evidently the differences between these transi- 
tion metal compounds are too large for a systematic relationship between ks and V, 
even though these also had the NaCl structure. Grimvall and Fernandez-Guillermet 
correlated the thermodynamic Debye temperature to the averaged mass and volume 
of a series of sp 3 -bonded compounds with the sphalerite (ZnS) structure [40]. The 
correlation included InSb, GaP, GaAs, and the elements Si, Ge, and Sn were added 
to the plot. They found a relationship of the form 0d [K] = 1500 a _1 M -1 / 2 [nm _1 
(at.wt.) -1 / 2 ] (where, for example 0d = 600 K corresponds to a~ l M~ 1 / 2 = 0.4 in 
these units). Adding volume differences to these types of correlations for different 
compounds having the same structure accommodates somewhat their differences in 
bonding and bond stiffness. Suggestions that there should be a correlation between 
interatomic forces and atomic volume, assuming sufficient similarity in chemistry 
and structure of the compounds, have been around for many years [41]. 

Quantitative analyses of the effects of bond length on bond stiffnesses were per- 
formed by Sluiter, et al. on Al-Li [42] and by van de Walle and Ceder on Pd-V [43]. 
They tested the idea of transferring interatomic forces from one material to another, 
with the goal of predictive power without major calculations or measurements on 
every material. There is some subtlety to the concept of transferrability of force 
constants because different crystal structures have different sets of allowed tensorial 
force constants. Force constants must be identically zero if symmetry operations 
of the crystal transform them into their negatives. Different symmetry operations 
for different crystal structures therefore make it difficult to transfer tensorial force 
constants from one structure to another. This is not a problem for axial force models, 
which involve a central-force "stretching" force and a "bending" force, and axial 
force models often give adequate descriptions of lattice dynamics. In the "bond- 
stiffness-versus-bond-length" approach, only first-nearest-neighbor (INN) bonds are 
transferred between different structures [9, 44]. It is also observed that only the 
stretching (longitudinal) forces between atoms have a significant dependence on the 
bond length. The bond bending forces are found to be smaller and less dependent 
on bond length. 

Van de Walle and Ceder performed density-functional theory calculations of bond 
stiffnesses for a set of fec-based Pd-V structures with different interatomic distances. 
Some results, presented in Fig. 6, show the longitudinal (stretching) interatomic 
forces between V-V, Pd-Pd, and Pd-V INN pairs of atoms. The points are from full 
calculations for individual Pd-V ordered compounds of different compositions and 
different lattice parameters. The points tend to fall rapidly with bond length, as 
expected when the interatomic forces soften with expansion. Approximately, there 
is a range of about 10% in the interatomic separation where the force constants 
vary from their maximum value to zero. Also approximately, the stretching stiffness 
has a linear relationship with bond length. The effects of bond length on bond 
stiffness of Fig. 6 are large, so a Pd-Pd bond cannot be expected to have the same 
stiffness in any local environment of a random a solid solution, for example. The 
bond-stiffness-versus-bond-length model is useful for semi-quantitative arguments, 
and when calibrated for specific elements by ab-initio calculations, for example, this 
model can be used for comparing phonon entropies of different alloy phases. The 
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bond-stiffness- versus-bond-length model would reduce to the bond proportion model 
if there were no variation of bond stiffness with interatomic distance. 

5.2. Extending the Bond Proportion Model 

Consider how a bond-stiffness- versus-bond-length relationship alters the free energy 
(Eq. 60) and the critical temperature (Eq. 62) of the bond proportion model for 
ordering transformations. The free energy F{L) of Eq. 60 is with reference to the 
disordered solid solution, and is zero for L = (i.e., F(0) = 0). Consider again the 
case of Bl ordering (NaCl prototype) of an AB alloy on a simple cubic lattice, as 
discussed following Eq. 62. In this structure, all INN bonds are between A-B pairs, 
and there is only one bond length. It is the disordered solid solution where the 
distortions of bond lengths occur. This is a complex problem owing to the different 
relaxations in the different local environments. We can get to the essence of the 
problem by using the same Griineisen parameter for all bonds (Eq. 65), and by 
simplifying the three curves of Fig. 6 by assuming that the bond stiffnesses are 
zero at the same bond lengths for A-A, B-B, and A-B pairs. In this case, all three 
frequencies {ujaa, <^bb> <^ab} change with volume as: 

"V = < f 1 " 1 S Z) . (67) 



v \ V° J 

where the superscript denotes a reference value for the ordered structure. For small 
5V/V 

ln(c^ln(4)- 7 ^. (68) 



To the free energy of Eq. 60 we add the term 

2 \T)AT ™ _ _ 

yo ' 



-TSoruen = k B T (1 - L 2 ) V N 7 ^ , (69) 
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which vanishes when L = 1, but grows with disorder as the number of R-W pairs 
increases as 1 — L 2 . The critical temperature of Eq. 62 becomes: 



Consider a typical case where the disordered phase has a lattice parameter 1% 
larger than the ordered phase (so 5V /V° ~ 0.03). Assuming 7 = 2, the second 
term in the denominator of Eq. 70 is +0.36. The critical temperature is suppressed 
substantially when the average bond length in the disordered solid solution is 1% 
longer than for the ordered phase. 14 Because we assumed the same Griineisen pa- 
rameter for all bonds, this effect of bond stiffness versus bond length is independent 
of the effect on T c from stiffness differences of the different pairs of atoms. In a more 
complete treatment these effects are inter-related. 

5.3. Chemical Effects on Bond Stiffness 

A closer examination of Fig. 6 shows that the bond stiffnesses for V-V and Pd-V 
pairs have considerable scatter about the linear fit. Some of this can be accounted 
for by a dependence of bond stiffness on chemical composition (discussed further 
in Sect. 20.1). For example, the addition of B atoms to a host crystal of element A 
will fill or deplete the electronic band structure of A, causing changes in the density 
of electron states at the Fermi level. The electron screening of the A-ion cores is 
thereby altered, changing the bond stiffnesses of A-A pairs. Such effects may scale 
with the electronegativity difference between A and B atoms, and features of the 
electronic structure such as the density of states at the Fermi level. These chemical 
effects are probably less significant for alloys of isoelectronic elements. 

Liu, et al., studied the stretching forces in fee aluminum alloyed with early tran- 
sition metals [45] , and presented a number of graphs such as Fig. 7. The stiffness of 
Al-Al bonds decreases with Zr concentration. This effect is big, and is comparable 
to the change of bond stiffness versus bond length at one alloy composition. The 
Pauling electronegativity of Al (1.61) exceeds that of Zr (1.33), so alloying should 
cause substantial changes in electronic structure. Liu, et al., argue that when Al is 
alloyed with Zr, electrons are transferred into the stronger Al-Zr bonds, stiffening 
them at the expense of Al-Al bonds. The stretching forces between Al-Zr and Zr-Zr 
pairs are much less affected, however. These effects of alloying on bond stiffness in 
Al-Zr are discussed further in Section 11. 



6.1. Harmonic Heat Capacity 

First consider a harmonic model where Cy is obtained with a fixed set of oscillator 
frequencies, or a fixed phonon DOS. The temperature derivative of the Planck 



k B T c 




(70) 



6. Heat Capacity 



1 In the less common case where the lattice parameter of the disordered alloy is smaller than the 
ordered, T c would be raised. 
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Figure 7. Nearest-neighbor stretching forces between Al-Al pairs in fee Al-Zr alloys with Zr 
concentrations as labeled. After [45]. 

occupancy distribution (Eq. 7), weighted by the energy of the mode, s^, shows that 
each mode of energy £j contributes to the heat capacity at T the amount 



In the high temperature (classical) limit where k-^T 3> £j, the heat capacity becomes 
1 /cb per mode. 15 Equation 71 is consistent with Eq. 6 if we integrate over T for all 
modes 



(Usually we integrate against the phonon DOS, g(e), to obtain SVvibC^)-) 

A Debye model predicts some general features of the heat capacity. Figure 8a 
shows heat capacity curves for two materials with similar Debye temperatures. The 
curves look similar at low T, and are essentially identical at very high T. Their 
most pronounced difference is found at approximately one-quarter of the Debye 
temperature. The integrand of Eq. 72 has a maximum below about a fifth of the 
Debye temperature, so experimental studies of vibrational entropy by heat capacity 
should include this range of temperatures. The difference between these two curves, 
shown in Fig. 8b, is typical of the shape of a differential heat capacity curve obtained 
from real phonon DOSs (see curve "DOS 11 K" in Fig. 26), although there are 
differences in detail. The Debye model assumes a phonon DOS that is a quadratic 
function of e, and increases monotonically up to a cutoff energy £d- The actual 
DOS often has a more abrupt rise at low energies that comes from the lower sound 
velocity of the low transverse modes, which are also weighted at their cutoff energy 

15 Even with phonon softening this could remain true. In the high-temperature limit, with phonon 
occupancies of feT '/e, per mode, the softening of the modes gives an increase in occupancy pro- 
portional to 1/Ei, but a decrease in energy proportional to Ei. The product of these two factors 
produces a heat capacity C'v at high temperatures that is unaffected by mode softening, so CV can 
remain as 3£;B/atom at high temperatures. This is not so if there are other phenomena at work, 
such as a buildup of elastic energy as a crystal expands. 




(71) 




(72) 
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by the effects of dispersion (e.g., from a high DOS at the Brillouin zone boundary). 
Plots of £d versus T typically show a dip and a rise from a mimimum at roughly 
a tenth of the Debye temperature #d = £dAb> but this varies from material to 
material. 

6.2. QUASIHARMONIC THERMODYNAMICS 

More typically measured by calorimetry is the amount of heat going into the solid 
at constant pressure, C p , defined as 

C P (T)=T^\ . (73) 

' p 

The difference between C p and Cy is a classical thermodynamic relationship 

C p -C v = 9Bva 2 T , (74) 

where B is bulk modulus, a is the linear coefficient of thermal expansion, v is specific 
volume, and 9 is the square of three dimensions. Equation 74 can be derived from 
formal thermodynamic manipulations, aided for example by techniques of Jacobian 
transformations. 

A physical derivation of Eq. 74 can be performed by considering the free energy 
of a crystal having only 1) phonons, and 2) thermal expansion. There are energy 
terms from the phonons and from the elastic expansion, and an entropy from the 
phonons (an average thermal expansion does not generate entropy itself). The free 
energy is 



F(T) = E elas (T) + E ph (T) - TS ph 



(75) 
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The elastic energy of thermal expansion, using SV/Vq = 3aT is 

E clas - ^B— , 



-^clas 



-BV a 2 T 2 , 



The phonon energy in the classical limit is 

E ph (T) = 3Nk B T . 
In the high-temperature limit, the phonon entropy is (cf., Eq. 15) 



(76) 
(77) 

(78) 



Using the same approach as Eq. 68 with Uj = u;°(l — j5V/Vq ) and SV/Vq = 3aT 
for the thermal expansion 



S ph (T) 



3N 



k B T 



+ 9Njak B T 



S° h (T) + 9N 7 ak B T = 5° h + S q , 



(79) 

(80) 
(81) 



where the first term, Sp h (T), is the harmonic phonon entropy with phonon frequen- 
cies characteristic of T = 0, and unchanged with temperature. Using Eqs. 76, 77 
and 81 in the free energy expression for the expanded crystal with phonons, Eq. 75 



F(T) = ^BV a 2 T 2 + 3Nk B T - TS° h (T) - 9Na~ 1 k B T 2 



(82) 



First we obtain the equilibrium thermal expansion coefficient, a, by calculating 



dF(T) 
da 



= 9BV aT 2 - 9N-/k B T 2 , 



a 



3BV 



(83) 



where we have used the classical result Cy = 3Nk B . 1(i If the electronic entropy 
depends on volume, an additional term, C c \y 7 e i/(3BVo) is added to the right-hand 
side of Eq. 83. 

Substituting Eq. 83 for the thermal expansion into 82 for the free energy 



F(T) = F°(T) - -BV a 2 T 2 , 



(84) 



16 The effect of phonon energy is small as discussed in Footnote 15, although there is also a small 
effect from a reduction in zero-point phonon energy as the solid expands. 
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Figure 9. Illustration of a relatively large effect of thermal expansion on the heat capacity, using 
models of a cerium alloy. Inset is the phonon DOS at K, from which the "harmonic" heat capacity 
curve was obtained (3fcB/atom is the high-temperature limit of the harmonic model). Using a large 
a and an experimental B in Eq. 83, a 7 was obtained from Eq. 83 that predicted the phonon DOS 
would soften by 6.5 % from to 300 K. The elastic energy from thermal expansion was added to 
the "harmonic" curve to obtain the "quasiharmonic" curve. 

where F°(T) is the harmonic free energy at finite temperature with zero thermal 
expansion (a = 0). We use Eq. 83 to express the quasiharmonic entropy S q of Eq. 
80 or 81 as 

S q = 9BV a 2 T (85) 

(consistent with differentiating Eq. 84 with respect to T). Finally, if we substitute 
Eq. 85 into Eq. 73, we can recover Eq. 74. Typical effects on heat capacity are 
shown in Fig. 9. This S q of Eq. 85 should always be considered when comparing the 
vibrational thermodynamics of different materials at elevated temperatures. 

6.3. Anharmonic Heat Capacity 

In mechanics, the word "anharmonic" describes any oscillator with generalized 
forces that deviate from linearity with generalized coordinates. In vibrational ther- 
modynamics, "anharmonic" is used more restrictively. The word "quasiharmonic" 
accounts for the effects of thermal expansion on phonons just considered (e.g., Eqs. 
74, 85), effects that are not strictly harmonic. An implicit assumption of quasi- 
harmonic theory is that an oscillation has a well-defined, long-lived, harmonic-like 
frequency that changes modestly and predictably with changes in T and P. 

Equation 83 is useful for separating anharmonic from quasiharmonic behavior. 
Suppose we know (by inelastic scattering measurements for example) the thermal 
softening of the phonon modes, from which we calculate an average Griineisen pa- 
rameter 7. With the phonon DOS we can also calculate Cy (or its molar quantity), 
and conventional measurements can provide B, Vq (or molar volume v), and a, 
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accounting for all unknowns in Eq. 83. If Eq. 83 proves to be true with the 7 
measured from the phonon softening, we say that the solid is "quasiharmonic." 

Anharmonic behavior is identified as phonon softening or stiffening that is incon- 
sistent with the Griineisen parameter 7 needed for equality in Eq. 83. Anharmonicity 
involves phenomena beyond those of independent phonons. Higher-order phonon- 
phonon interactions are the source of anharmonic behavior as discussed in Sect. 
18.2. Section 20 discusses electron-phonon interactions, which also vitiate Eq. 83. 17 
Substantial deviations from quasiharmonic behavior are known for Cr [47], V [48], 
Ce [49], Mo, W, Si, Ge [50], and anharmonic or non-quasiharmonic effects could be 
important in many other systems. 

6.4. Thermodynamic Entropy 

The temperature dependence of the entropy 



contains harmonic, quasiharmonic, and anharmonic contributions, along with non- 
harmonic effects from electron-phonon interactions. For the harmonic part of the 
vibrational entropy we use Eq. 14, with a phonon DOS, go(e), from low temperature 
measurement or calculation (ideally T = OK). The harmonic part contributes to 
the first term on the right-hand-side of Eq. 86, i.e., dS/dT)y. Harmonic phonons 
undergo no change in frequency with T or V. The temperature dependence of Sh(T) 
originates with the Planck occupancy factors n(e,T). 18 

Quasiharmonic phonons have frequencies that depend on volume only, and make 
a contribution through the second term on the right-hand-side of Eq. 86. At a fixed 
volume, however, they behave as harmonic oscillators. Their frequencies can change 
with temperature, but only because thermal expansion alters the volume of the solid. 

Anharmonic entropy, S a , may contribute to either of the two terms on the right- 
hand side of Eq. 86, but especially the first, a pure temperature-dependence of 
the vibrational entropy at fixed volume caused by changes in the interatomic force 
constants with temperature. Consider for example a quartic term in the interatomic 
potential, which alters the force constants with larger thermal displacements. The 
quartic term does not affect the thermal expansion, but it does change the phonon 
frequencies. It causes a type of "phonon-phonon interaction," so called because 
the existence of phonons affects the energy to create new phonons, as discussed 
in Section 19. Other temperature-dependences of force constants at constant vol- 
ume may arise if the nature of the bonding changes as atoms are displaced further 
from their equilibrium positions, or perhaps thermal excitations of electrons alter 
the sensitivity of electronic energy to atom displacements (Sections 17, 18.3 and 

17 It remains a challenge to differentiate between these sources of non-harmonic behavior [46], 
so sometimes the effects of electron-phonon interactions are called "anharmonic." More often 
"anharmonic" is reserved for phonon-phonon interactions only. 

18 For harmonic solids at high temperatures, differences in the phonon DOS will cause a difference 
in entropy for different phases, but this difference in entropy does not change with temperature. 
(At high T, S increases with T by the same amount for each phase.) 




(86) 
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20 describe such effects from adiabatic electron-phonon interactions). The second 
term in Eq. 86 is associated with changes of the interatomic force constants with 
volume. The quasiharmonic contribution is expected to account for most of it, but 
the practice of using harmonic entropy expressions, while correct to first order [46], 
is not fully verified. The cubic term associated with thermal expansion also shortens 
phonon lifetimes, causing energy spreads of phonon linewidths. This may not alter 
the vibrational entropy to first order, but again the robustness of the harmonic 
approximation is not fully tested. 

Analogs exist for magnetic dependencies, as can originate from a coupling be- 
tween magnetic and vibrational degrees of freedom [51]. Invar behavior is an ex- 
perimental example, but perhaps this is best understood as a volume-dependent 
electronic transition rather than a magnetic one [52, 53, 54]. 

7. Calorimeters 

The word "calorimeter" covers a broad scope of instruments for measuring large 
or small transfers of heat to samples at low or high temperatures. From these 
measurements heat capacity can be determined, or heats of phase transformations 
or chemical reactions in the sample can be found. 19 Phase transformations and 
chemical reactions may involve changes in configurational entropy, however, so for 
isolating vibrational entropy it is usually better to measure heat capacity data at 
lower temperatures where no changes in local atomic arrangements are underway. 
(An exception is the martensitic transformation discussed in Sect. 15.) 

Adiabatic calorimeters have performed heat capacity measurements for a century, 
dating from Nernst and Eucken [55, 56, 57, 58, 59, 60]. The sample and a heater are 
thermally isolated from their surroundings by vacuum and an adiabatic shield that 
is maintained at approximately the same temperature as the sample. The latter 
minimizes radiative heat transfer, and a thin mechanical suspension such as silk 
thread minimizes heat conduction. Measured electrical energy is provided to the 
sample through a heater, and the sample temperature is measured by thermometry. 
After equilibration the temperature rise is measured, from which the heat capacity 
of the sample and its holder are obtained. Ideally, the sample weight should be large 
(lem 3 is relatively small), minimizing background corrections. Heat capacities can 
be measured to better than 1% accuracy from typically 10-300 K with adiabatic 
calorimeters, but there are many instrument designs and operating ranges. A varia- 
tion on this method is drift calorimetry. The sample is connected to a thermal bath 
through a weak thermal link, causing its temperature to drift. Although the rate 
of drift should be proportional to the heat capacity of the sample, more accurate 
data are usually obtained by applying pulses of electrical energy, and identifying 
the extra temperature change over the baseline drift. Another approach is used for 
heat pulse, or relaxation calorimeters, which also have thermal links between the 
sample assembly and a thermal anchor. After a pulse of electrical energy raises 
the temperature of the sample, its temperature decays towards that of the thermal 

19 Proposals appear with some regularity to use accurate heat capacity data to determine the 
phonon DOS. These prove to be too optimistic, although some lower moments of the phonon DOS 
can be determined this way. 
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anchor. The rate of decay (sometimes modeled as an exponential function or a sum 
of two exponentials) is used to obtain the heat capacity of the sample, which can 
be proportional to the characteristic time of the decay [61, 62, 63, 64]. A related 
technique, which also has advantages for small samples, is to apply an oscillating 
heating signal to the sample, and measuring the thermal response at the frequency 
of the heating signal using a lockin amplifier to minimize noise [65]. 

Differential measurements can be highly accurate, and are sometimes useful for 
measuring differences in vibrational entropy of different materials. A differential 
scanning calorimeter heats two samples, perhaps matched in their composition and 
mass, but with different crystal structures. The electrical heating is performed with 
feedback control so that the temperatures of both samples increase at the same 
rate. Differences in heat capacity are measured by the differences in electrical power 
required to maintain the same heating rates. This method of differential scanning 
calorimetry has, in principle, the advantage that it nulls the background from the 
sample holder, but background issues prove difficult when highly accurate measure- 
ments are required. The method also has limitations in the range of temperatures 
that can be achieved, because if there is no adiabatic isolation, the electronic heating 
system needs to provide considerably higher power at higher temperatures. A simpler 
differential measurement technique, differential thermal analysis, places two samples 
in a common thermal environment, such as side-by-side in a furnace. The furnace 
temperature is increased uniformly, and the temperature difference between the two 
samples is measured to obtain their differential heat capacity. These differential 
measurements can be accurate to 2% of the total heat capacity, but like all methods 
of calorimetry, it is a challenge to obtain high-quality results over a wide range of 
temperatures. 



8. Neutron Scattering 

8.1. Elastic and Inelastic Scattering 

Figure 10 shows the geometry, wavevectors and position vectors for neutron scat- 
tering. The plane wave incident from the left, *i n (r , ',t), has the form 

* in (r^) = e^-"'*) . (87) 

For scattering from one nucleus (which is essentially a point compared to typical 
neutron wavelengths), a spherical wave, * sc (r, t), emerges from r", and reaches the 
detector at r at time t' > t 

_1 e i[fc f |r-f |-u; f (f-t)] 
%c(r, t') = — — . (88) 

4vr |r — r | 

The difference, f—r 1 , is the distance from the nucleus to the detector. Because Eq. 
88 gives the contribution to the wavefunction at position f caused by a 5-function 
at r*, with its prefactor it is the Green's function of the Schrodinger equation of the 
neutron scattering problem [34]. 

In the first Born approximation [34], the amplitude of the scattered wave at 
the detector, * sc (r, t'), is proportional to the product of: 1) the amplitude of the 
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kj detector 




Figure 10. Wavevectors and position vectors for neutron scattering. 



incident wave, ^ m (Eq. 87), 2) the Green's function of Eq. 88, and 3) the potential 
that scatters the neutron out of its initial state 

m r - i[fcf|r-f|-a;f(t'-t)] 

*"< ? > <> " S / W "*"""' |r- r -| dV ' < 89 > 

For nuclear scattering, relevant to phonon measurements, we use the "Fermi pseu- 
dopotential," V nuc (r f ) = Airh 2 / (2m n ) 6(5(r y ) which places all the potential at the 
point nucleus, where b is a fixed length. 20 Fermi pseudopotentials are placed at 
crystal translation vectors, {Rj}, accounting for all TV nuclei in the crystal 

h 2 N 

^) = 4 ^EW-*i)- (9°) 

This sum of (^-functions transforms the integration of Eq. 89 into a sum over points 
at Rj . For laboratory distances (large compared to the size of the scatterer) ,r^>Rj, 
the outgoing spherical wave can be treated as a plane wave, allowing rearrangement 

„i(fcf -r— oj{t') _ _ - 

q> sc (f, = - — h i e i(fci " fef) -^ e i(wf - Wi)t . (91) 

\ r \ j 

The factor ^(kfr-uift') gj ves ^ e precise phase at the time and place of neutron 
detection. It has modulus unity, and we ignore it because we measure the intensity 
of scattering. The relative phases of waves from all nuclei are essential, however, 
because they account for interferences. 

We define the change in neutron wavevector, Q = k\ — k{ (and the momentum 
transfer from the neutron to the sample is HQ). Using Eq. 18 for Rj, and allowing 
the scattered wave to emanate from the nucleus at later times 

i£ sc (f) = I b K e iQ-(n+r«) e ifci-« iK (0) e -ifcr« iK W e i(^-^i)t ; ( 9 2) 



r 



where we note that because all unit cells are identical, b depends on the basis vector 
only. It is easiest to understand the case when the scattering is a fast process, such as 



20 Typically of a nuclear scale of 10 12 cm or so, sometimes with an imaginary component to 
account for absorption. 
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for x-ray scattering, where an x-ray traverses the crystal more rapidly than typical 
vibrational periods of atoms. Incoherent neutron scattering is also fast, where the 
thermal neutron is scattered from individual atomic nuclei without concern for phase 
relationships between nuclei. This motivates a first approximation; ignore the time 
dependence of the potential by setting ui K (t) = ui K (0), and take the time average 

T 

Vf c {Q,r) = 1^ 6 K e i<5 '^ + ^ e i( ^ (0) i j e^'^dt , (93) 

l,K 

*s'(Q, ^ = \ E b - e^i+r.) e iQ-« iK (o) §{uJ{ _ Wi) ^ (Q4) 

l,K 

The 5-function enforces no change of frequency or energy, so Eq. 94 pertains to 
elastic scattering, where the modes of atom vibrations are unchanged. We need to 
perform an ensemble average over the different values of ui K (0), which evaluates as 21 

Equation 94 for elastic scattering becomes 

P -W _ . P -W 

ni(Q, ?) = — E h - eiQ ' rK E eiQ ' n = — hq) s(Q) . (96) 

r K I r 

Here we separated the sums over basis and lattice vectors into a structure factor, 
F{Q), and a shape factor, S{Q). Equation 96 can be used to calculate the neutron 
diffraction pattern from a crystal, giving the average static structure. The factor e~ w 
is the "Debye-Waller factor," which becomes appreciably less than 1 when the atom 
displacement u K is comparable to 1/Q ~ A/27T, the wavelength associated with the 
scattering angle and energy. This is consistent with the origin of the Debye-Waller 
factor from the destructive interference of scattered wavelets that are emitted from 
different spatial positions of the displaced nuclei. Note that the ensemble average of 
Eq. 95 does not reveal the energy of an atom vibration, only how displaced atom 
positions alter the phase. 

Next consider a process of scattering that is slow compared to the atomic vi- 
brational period. Such is the case for coherent scattering of neutrons by phonons, 
and for nuclear resonant scattering of photons by phonons. For clarity, we consider 
the scattering by one phonon of wavevector k. Equation 28 shows that u t = 
{t ;K ^(0)e~ lw fc*. We need to take a time average of the displacement factors in Eq. 92 

(e^KfcW e _i ^ f '"i«fc ( * ) ) t = (e^'^-fc^ e"^*""!^ ^ 1 ^*^ . (97) 

Coherent phonon scattering involves correlations of displacements u [K ^(t) at dif- 
ferent atom positions at different times. The intensity of the coherent inelastic 

21 An approximate demonstration of this is by expanding exp(iQ ■ u) — 1 + iQ ■ u + ... and taking 
either an ensemble average or a time average. The average of the second term, iQ ■ u is zero owing 
to positive and negative excursions of u, leaving the quadratic term, consistent with Eq. 95. The 
average of the subsequent cubic term also vanishes, conveniently, but unfortunately the quartic 
term has the wrong sign. Nevertheless, it can be proved that Eq. 95 is exact if the probability 
distribution of u is Gaussian about the origin, and this is often reasonable for lattice dynamics. 
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scattering proves to be the Fourier transform of the Van Hove space-time correlation 
function [66], and is a rich approach to studying dispersive excitations in solids. 
The proper quantum mechanical treatment of time and temperature averages is 
beyond the scope of this review, but Eq. 98 below is the proper result for incoherent 
scattering. Furthermore, it is possible to average the coherent scattering over all Q to 
recover this incoherent result [67, 68, 69], and the incoherent approximation serves us 
well in the following analysis of multiphonon scattering. The loss of phase informa- 
tion by averaging the coherent scattering over all Q brings us back to consideration 
of a fast scattering process. 

8.2. One-Phonon and Multiphonon Scattering 

Equation 96 is the elastic scattering only. It does not conserve the nuclear scattering 
cross section. Nuclear interactions do not necessarily depend on the displacements 
of atoms, so we need a compensating factor e + ^' M ) ^ to conserve the cross-section 
and account for the rest of the scattering intensity. Expanding this new factor, we 
calculate the intensity as ^* c ^ sc , and for clarity we normalize the intensity 

<*s*c*sc>iV = ^X>- <( ^ )2) 

x (l + ((Q-u K ) 2 )+ 1 -(((Q-u K ) 2 )) 2 + ...) (98) 

The first and second terms in parentheses of Eq. 98 are the fractions of elastic (0- 
phonon) and 1-phonon scattering. The third term with the \{{{Q • u K ) 2 )) 2 accounts 
for 2-phonon scattering processes. Only the 1-phonon scattering proves useful for 
obtaining a phonon DOS, but it is important to have an understanding of the 
higher-order terms (called "multiphonon scattering") for designing an experiment, 
or performing quantitative corrections of experimental data. 

A more detailed treatment of coherent 1-phonon scattering includes the ther- 
mal factors for phonon populations, giving the double-differential scattering cross 
section 22 

d 2 a kt (2tt) 3 ^ n(k) + 
dQdE ~ h 2V ~t 

9,k k 

x|V-^= e-«^) 2 > [Q ■ e K (k)]e 1 ^ * 

x S(Q + k ± g) 5{uj - loj:) , (99) 

where the phonon occupancy factor, n(k), is the Planck distribution of Eq. 7 at the 
temperature of interest and ei = huur. The signs in the factor n(k) + \ ± ^ in Eq. 99 
are such that it is n(k) for phonon annihilation, and n{k) + l for phonon creation - it 

22 A new factor is the ratio fcf/fe, which expresses the effect on flux caused by the rate at which 
neutrons leave the sample. Compared to an elastically-scattered neutron, fewer neutrons per second 
will pass into an increment of solid angle dQ if they are slowed to smaller values of k{. Also, the 
use of solid angle eliminates the need for the factor r~ 2 . 
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is always possible to create a phonon, even at T = when no phonon excitations are 
present. The 5(Q + k±g) and 5(uj — uj(k)) enforce the conservation of momentum and 
energy on the neutron scattering process. The phonons on different branches must 
be considered separately in Eq. 99, and it is traditional to add a "branch index," 
sometimes denoted j, to u(k), since more than one ujj may correspond to a specific 
k (cf., Eq. 29). 

Why should the expansion of an evidently classical exponential in Eq. 98 produce 
the series of neutron-phonon interactions? In the classical limit, the scattering of 
every neutron is inelastic, generating heat composed of a number of phonons with 
total energy, E, the energy transfer from the neutron. For classical scattering, an 
atom of mass m will recoil with velocity v = HQ/m, with momentum transfer HQ, 
and energy transfer 

ti 2 O 2 

Bt-OJL. (100) 

In the quantum limit where Q is small, the energy recoil may or may not occur, 
depending on whether phonons are created or not. The total scattering, S(Q,E), 
is the sum of components, J2^=o Sn(Q, E), from neutrons scattered after creating 
different numbers, n, of phonons in the sample. 23 The thermal energy is 

k B T = 2muj 2 (u 2 ) (101) 

for one harmonic mode of frequency to. From Eqs. 100 and 101 

<?V> - § g U02) 

This is the quantity in which the multiphonon expansion is performed (Eq. 98). 
From Eq. 102 we see that the relative intensities of the sequence of phonon scattering 
terms scale with 1) the ratio of the recoil energy to the oscillator energy, and 2) the 
ratio of temperature to the oscillator energy. When the first ratio is small, quantum 
mechanics requires that some recoils excite the oscillator, but most do not. The 
same is true for temperature (and at low temperatures we of course have low heat 
capacity because the oscillator is usually not excited). 

It remains to get the spectral shape of each order of the multiphonon scattering. 
The spectrum for one-phonon scattering weights more heavily the low-energy modes 
because they have larger amplitudes of motion, providing a factor of g{e)/e. The 
number of phonons is the Planck distribution l/[exp(e/k B T) — 1], and the two factors 
provide the shape of the 1-phonon profile, Ai(e): 

When two phonons are created simultaneously, the total spectrum of energies is 
the convolution of the 1-phonon profile with the 1-phonon profile. 24 The n-phonon 



23 Phonon annihilation is handled by extending the range of E to negative values for each 
S n {Q,E). 

24 Consider each phonon excitation to be a random variable with a probability distribution of 
Ai. The sum of two random variables has a distribution that is the convolution of the probability 
distributions Ai * Ai . 
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profile is the convolution of another 1-phonon profile with the (n — l)-phonon profile 

oo 

A n (e) = A 1 * A n ^ = J Ai(e - e') A n ^(e') da' . (104) 

— oo 

Starting with A±, we can generate the spectral shapes of all orders of multiphonon 
scattering by the systematic application of Eq. 104. 25 The total scattering is the 
sum of these spectral profiles, weighted by the corresponding terms of Eq. 98 

00 (2W) n 
S(Q, E) = J2 exp(-2W)^-±- A n (E) , (105) 

77-. 

where we have defined S(Q,E) as the double-differential cross-section d 2 a/dQdE 
without prefactors. 26 

In much experimental work, we seek to isolate the single scattering profile, A±(E), 
from the measured S(Q,E), and then determine the phonon DOS g{E) by Eq. 103. 
One way to correct for multiphonon scattering is by iteration on an initial guess of 
the phonon DOS. With an initial guess at the DOS or Ai(E) (often the S(Q, E) with 
a constant background removed), the various ^(-E 1 ), A${E), A^(E)... are calculated 
and subtracted from S(Q,E) to give a better A±(E). Iteration is complete when 
there is a good match between the measured S(Q,E) and the S(Q,E) recalculated 
with Eq. 105. When the multiphonon scattering is weak, perhaps only two iterations 
are necessary. When the multiphonon scattering is strong, the procedure becomes 
difficult and slow. If there is strong multiphonon scattering, however, the experi- 
ment is moving into the limit of classical scattering, so isolating the single phonon 
spectrum becomes more challenging. 



8.3. Neutron Weighting 



Strictly speaking, the phonon DOS from the procedure outlined above is not the 
true phonon DOS, but rather the neutron-weighted DOS. For an elemental scatterer, 
the neutron-weighted DOS is identical to the phonon DOS. This is not the usual 
case for a sample that contains more than one type of atom. The neutron-weighted 
phonon DOS is 

9Nw(.e) oc Y<9d(e)exp(-2W d )eM2W)^±- (106) 
d md 

where the ^rf(e) are the phonon partial DOS functions of Eq. 34. The Debye- Waller 
factor is an explicit function of gd{z)- The term exp(2VF) is the average Debye- Waller 

25 A subtlety arises for coherent scattering from crystals that are sufficiently anharmonic so 
that the phonon lifetime is shorter than the characteristic time for neutron scattering. When, for 
example, one phonon decays into two during the scattering, there can be constructive or destructive 
interference between the neutron wavepackets scattered by a one-phonon process and a two-phonon 
process. This interference can alter the spectral shape. Although recognized some time ago [70], 
identifying this effect has been a challenge [71, 72, 73]. 

26 This S(Q,E), called the "scattering law," is normalized so that it is the double Fourier 
transform of the Van Hove correlation function G(f,t). 
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correction; it is calculated from the self-consistent neutron-weighted DOS. At low- 
temperatures, the factor exp[2(VF — W^)] is approximately unity. 

For the case where ad/rrid is the same for all species d, <?nw( £ ) ~ 9{ £ )> but in 
general obtaining the true phonon DOS from the neutron-weighted phonon DOS 
requires a full analysis of the lattice dynamics. This can be performed by simulational 
procedures or by some types of analytical methods. These are important analyses, 
but are mostly beyond the scope of this review. One example is presented in Section 
13.2. 

8.4. Direct-Geometry Fermi Chopper Spectrometer 

Chopper spectrometers have evolved with advances in electronics [74, 75, 76, 77, 78]. 
A Fermi chopper is a fast shutter located a fixed distance from the neutron source. By 
aligning a slot at the right time after a neutron burst leaves the moderator source, it 
selects neutrons having a particular velocity, and hence a particular neutron energy 
E = T^mv 2 . Its operation is best illustrated with the distance versus time diagram 
of Fig. 11. The slopes of the lines are neutron velocities. Neutrons of many velocities 
are emitted from the moderator. The figure depicts two neutron pulses from the 
moderator, separated in time by perhaps 1/60 sec. The Fermi chopper selects a 
narrow range of velocities, corresponding to a time window of tens of microseconds. 
Of course the narrower the window, the more neutrons are blocked by the chopper, 
and the lower the flux of neutrons on the sample. On the other hand, the narrower 
the window, the more precisely selected are the velocities, and hence better energy 
resolution is achieved. The Fermi chopper therefore controls the incident energy, the 
intensity, and the energy resolution of the neutrons incident on the sample. 27 

The sample usually transmits most of the incident neutrons without energy 
change, and most of the scattered neutrons are scattered elastically. There is a 
strong tendency for the lines in Fig. 11 to have unchanged slopes through the 
sample. Neutrons that are scattered inelastically, however, have lines with kinks 
at the sample, and have slopes after the sample that are either steeper or shallower 
than for the incident beam. 

A Fermi chopper spectrometer works by timing. The time of detection of the 
elastic neutrons is known from the velocity selected by the Fermi chopper, and 
the elastically-scattered neutrons give an intense peak in the time spectra from the 
detectors. The neutrons that arrive earlier have gained energy from the sample, and 
those that arrive later have lost energy to the sample. The neutron energy spectrum 
is obtained from the histogram of neutron arrival times. From the energy transfer 
and the scattering angle, the momentum transfer can be determined. 

8.5. Triple-Axis Spectrometer 

The triple-axis spectrometer, invented by Brockhouse [79, 80, 81], uses a steady beam 
of neutrons from the core of a nuclear reactor. After the neutrons are moderated to 
thermal energies, a monochromator crystal selects neutrons of a particular energy 

27 The moderator also affects energy resolution. Figure 11 shows emission of neutrons at two 
points in time, but there is actually a spread of emission times, especially for neutrons of lower 
energy. 
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Figure 11. Distance-time diagram for inelastic scattering by a direct geometry chopper spectrom- 
eter. Positions of the moderator, chopper, sample, and detectors are marked on the vertical axis. 
Two moderator pulses occur at points on the z-axis. 

and sends them in a beam to the sample. The sample is typically a single crystal, 
positioned with different orientations with respect to the incident beam. An analyzer 
crystal is used to select an energy of the neutrons leaving the sample, and directs the 
final beam into a detector. There are a number of modes of operation of a triple-axis 
instrument, but for measuring phonon dispersions along crystallographic directions 
it is typical to operate the instrument with a fixed momentum transfer, Q. The final 
energy can be selected by fixing the angle of the analyzer crystal with respect to 
the sample and detector. The energy transfer is varied by changing the orientation 
of the monochromator crystal, and moving the sample to remain in the incident 
beam. Maintaining a fixed momentum transfer requires changing the orientation of 
the analyzer system around the sample. A computer system is used to control the 
scan over E at constant Q. 

The triple-axis spectrometer selects one value of E and one value of Q at any one 
time. It is the instrument of choice for measurements on single crystals, where it 
is important to measure dispersions of excitations along crystallographic directions. 
Single-crystal measurements of phonon dispersion curves can be interpreted directly 
when phonon peaks are evident above the background. Fitting the dispersion curves 
to Born-von Karman models has been an important method for determining inter- 
atomic force constants. With these force constants it is then possible to calculate 
the energies of phonons of other k, and obtain a phonon DOS. 
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9. Ab-initio Methods 



9.1. Density Functional Theory 

A large branch of computational materials science owes its predictive power to 
developments in density functional theory (DFT), in conjunction with advances 
in computing power. There are excellent reviews and books on density functional 
theory and its applications to materials, some of which offer a broader context of 
condensed matter theory [82, 83, 84, 85]. This rich field is far too extensive to 
cover here. The techniques for calculating phonons from first principles also deserve 
a much more extensive presentation [9, 86, 87, 88, 89, 90, 91, 92, 93, 94]. This 
brief overview of some key methods and concepts is presented with apologies that 
numerous intellectual achievements in this rapidly-moving field cannot be covered. 

Density functional theory elevates the importance of the electron density with 
respect to the electron wavefunctions, relying on a theorem that all properties of an 
electronic system can be obtained from the ground state electron density [95, 96]. 
A functional relationship must exist between the density and the physical quantity 
of interest, but finding this relationship is a challenge. The successful approach 
developed by Kohn and Sham [97] assumes that the ground state electron density of 
a system of interacting electrons can be equal to that of a non-interacting system of 
electrons, with an interacting electron density. Working with independent electrons 
makes calculations practical, freeing us from time-ordered-fermion field operators 
for electron-electron interactions. In particular, the electron exchange and Coulomb 
correlation energy, which are nonlocal multibody interactions, can be successfully 
approximated with functionals of the local density, or the density and its local 
gradients, for example. The density is obtained from the electron wavefunctions, 
so computational algorithms require iteration to self consistency of the potential 
and the wavefunctions. The local density functional theory has been surprisingly 
successful, although problems occur for strongly-correlated electron systems. 

9.2. First-Principles Phonon Calculations 

Principles for calculating phonon dispersions from the electronic energy were out- 
lined in the classic work of Ziman [98] and others [99]. Two approaches are used. The 
direct methods impose displacements of atoms, and calculate energies of distorted 
configurations. In the supercell method, pairs of atoms are displaced by u and u' 
along directions of symmetry, for example, and from the energies of the distorted 
structures, E(u,u'), the force constants 3>/ K £' K / of Eq. 27 can be calculated as 



(cf., Eq. 21). With these ab-initio force constants, the lattice dynamics can be cal- 
culated [88] . The Hellman-Feynman theorem shows that errors in the wavefunctions 
that cause only second-order errors in energy are responsible for first-order errors in 
forces, so energy convergence is stringent [9, 92]. With the frozen phonon method, the 
potential energy of an individual phonon is calculated directly by generating periodic 




(107) 
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distortions of a cell or supercell [100, 101]. Calculations with direct methods are 
versatile, and can account for unusual structures. Frozen phonons of different ampli- 
tudes are a natural way to account for anharmonic potential energies [102, 103, 104], 
although a good sampling of all phonons is challenging. When periodic boundary 
conditions are used, the cell should be comfortably larger than the longest-range 
interatomic force, and the phonon wavevector should be commensurate with the 
structure. With many atoms and electrons in a large supercell, however, the calcu- 
lations can become unwieldy, especially since the computation time often scales as 
the cube of the number of electrons. 

The second approach calculates the sensitivities to small density modulation, 
usually with a linear response theory [86, 87, 89, 105]. These are second derivatives 
of the energy, and the simplest approach to the problem is with second order per- 
turbation theory. For periodic crystals with both electrons and phonons expressed 
as Bloch functions, the density susceptibility function for the Fourier coefficient (k) 
of the electron density depends on the Fourier coefficient of the phonon modulation 
(k') as 

x&t - = 2EE (rt>*^), (108) 

where the intermediate electronic states \j) are unoccupied, and the states \i) are 
occupied. A considerable simplification is that k and k! must be related by a re- 
ciprocal lattice vector, reducing the number of terms in Eq. 108, although for more 
realistic interacting systems, the situation is not so simple. A variational approach 
for obtaining the energy to second order has been developed to obtain a linear 
response to periodic perturbations [89]. In this method, the dynamical matrix can 
be calculated directly [90]. It is also possible to recover real-space force constants 
by Fourier transformation of the dynamical matrix for appropriate selections of k. 

9.3. Molecular Dynamics 

Molecular dynamics simulations integrate equations of motion to follow in time 
the positions and velocities of all atoms. Many molecular dynamics simulations 
are classical or semi-classical, using forces obtained as gradients of model potential 
energy functions, for example. The atoms are allowed to move under Newton's law 
for a short time, and the forces are recalculated for the new configuration of atoms. 
These simulations can be extended beyond a microcanonical ensemble by allowing 
contact of a small simulation system with a reservoir, accounting for fluctuations of 
temperature or pressure, for example. Realistic thermostatic or barostatic control 
can be achieved by extending the dynamical system with degrees of freedom that 
represent the external system [106, 107]. Molecular dynamics is not restricted to 
harmonic forces, and large-amplitude displacements allow investigations of anhar- 
monic behavior. Accurate force fields that account for all detail of the local atomic 
environment can be a challenge. Such force fields can be obtained from quantum 
mechanics, but often they are empirical. 
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The velocity- velocity autocorrelation function (v(t)\v(0)) can be calculated from 
the trajectory of each atom in a molecular dynamics simulation. 28 The Fourier 
transform of the velocity-velocity autocorrelation function gives the phonon DOS 



In another approach, it is possible to obtain differences in thermodynamic quantities 
such as phonon free energy by parameterizing a Hamiltonian to allow the system 
to switch from one structure to another, and integrating along the path of this 
switching parameter. 

It is possible to solve for the forces on the moving nuclei by performing a full 
DFT calculation at each timestep of the simulation, but this is cumbersome. The 
Car-Parrinello algorithm makes practical quantum molecular dynamics simulations 
of greater complexity [109, 110]. It uses the electron degrees of freedom as dynamical 
variables in a set of coupled equations of motion for the nuclei and the electrons. The 
nuclei have a classical kinetic energy of M/2 [dui K /dt] 2 . The trick to the algorithm 
is that the electrons are given a fictitious kinetic energy of /x/2 f[dip/dt] 2 dr, with a 
mass \i that is an adjustable parameter (that helps with the stability of the numerical 
calculations). With this coupling of nuclear and electronic motions, constrained 
by the orthogonality of the electron wavefunctions, the electrons stay near their 
ground states during the simulation. The electron wavefunctions do not need to be 
recalculated at each step of the simulation, offering an improvement in efficiency. 
The Car-Parrinello algorithm is an adiabatic approximation. 



[108] 




(109) 



28 The notation indicates that the velocity history of each atom is multiplied by the velocity 
history shifted by t, and integrated over all time. 
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II. Vibrational Entropy of Materials at Lower 

Temperatures 



10. Earlier Investigations 



10.1. Polymorphism 

Polymorphism is the transformation of a material into a new crystal structure by 
changing T or P, without change in composition. The classic texts by J.C. Slater 
Introduction to Chemical Physics (1939) [111] and by F. Seitz Modern Theory of 
Solids (1940) [112] have explanations of the role of vibrational entropy in polymor- 
phic phase transformations of pure elements. An early experiment was reported by 
F. Lange in 1924, who measured the heat capacities of white and gray tin from 
nearly liquid helium temperature to ambient [113]. The difference in entropy was 
obtained with the Nernst postulate and by integrating the difference in heat capacity 



The enthalpy of the transformation was also obtained, and the critical temperature 
was predicted accurately from the vibrational entropy and the enthalpy [112]. It 
appears that the entropy of the tin transformation is mostly vibrational. It is big 
- white tin has a vibrational entropy 0.9fcB/atom greater than gray tin at room 
temperature. The roles of vibrational entropy in the polymorphism of sulfur, carbon, 
and iron are discussed semiquantitatively in Seitz's book [112]. 

Many elemental metals transform from close-packed structures (fee or hep) at 
low temperatures to the bec structure at high temperatures (and ambient pressure). 
These include metals with very different electronic structures, such as alkali met- 
als and transition metals. Zener [114] hypothesized that soft transverse phonons, 
intrinsic to the more open bec structure, would contribute enough vibrational en- 
tropy to favor the bec structure at high temperatures. Many years later, this idea 
was subjected to some analysis by Grimvall and Ebbsjo [115], using results from 
pseudopotential calculations [116, 117] and results from force-field calculations of 
lattice dynamics [118]. Their assessment showed that phonons in bec crystals were 
typically a few percent softer than in fee or hep, especially for the low transverse 
modes, supporting Zener's hypothesis. They further argued that this trend could be 
understood without anharmonic effects, although perhaps not quantitatively. 29 On 
the other hand, the fee and hep close packed structures seem to have very similar 
phonon energies, as shown by neutron scattering measurements on Co [119], for 
example. Likewise, the fee-based LI2 structure of C03V has a vibrational entropy of 
only 0.07 ks /atom larger than the ordered hexagonal hP24 structure of C03V near 



(cf, Eq. 72) 




(110) 







29 Their calculations, which assumed a constant volume per atom in both phases [115], arc 
approximate. 
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the transition temperature of 1320 K [120, 121]. On the other hand, for pure 57 Fe, 
nuclear resonant inelastic x-ray scattering (NRIXS) measurements of the phonon 
DOS at elevated pressure and temperature show only a small difference for the bec 
a-phase and the hep e-phase [122, 123]. Curiously, another NRIXS study has shown 
a large difference between the phonon DOS of the hep e-phase and the fee 7-phase 
[124], contrary to the trend for other metals. The phonon DOS of 7-Fe was softer 
than the DOS for e-Fe, giving the fee 7-phase the larger vibrational entropy and 
lower sound velocities. 

Electronic entropy also plays an important role in the hep-bec transformation 
[17]. Especially today, ab-initio methods are well suited for comparing vibrational 
free energies of polymorphic transformations of many elemental metals and com- 
pounds (e.g., [125]). 

10.2. Solubilities and Defect Concentrations 
10.2.1. Diffusion 

Vibrational entropy originates with atom displacements around equilibrium, which 
are small compared to the atomic displacements of diffusional jumps. When the 
times between diffusional jumps become very short, however, the distinction between 
vibrational and configurational degrees of freedom becomes blurred. This is certainly 
an issue when jump times are comparable to vibrational periods, or even if jump 
times are comparable to phonon lifetimes. This blurring of vibrational and diffusional 
phenomena occurs in liquids, for example, where some low-energy modes of vibration 
become unstable and atom movements occur. Likewise, soft-mode transformations 
relate vibrational degrees of freedom to structural ones as discussed in Sect. 15. 

Even in the common situation where the energetics of the atomic displacements 
are of order k-gT ~ 0.1 eV and the activation barriers for diffusion are of order 1 eV, 
there are relationships between vibrational entropy and diffusional processes. The 
activation free energy for diffusion, AG*, is related to the jump frequency, T, as: 

r(T) = r e - AG */ fc BT = ro e+ A S */k Be -AH*/k B T ? 

where Tq is an attempt frequency related to atom vibrations, and AS* is analyzed 
in terms of vibrations. Zener estimated AS* by how elastic constants soften with 
temperature [126]. 

The relationship between atomic vibrations and the entropy of diffusion, AS*, 
was developed formally by G. Vineyard, who assigned the origin of AS* to the 
vibrational entropy of all harmonic normal modes of the equilibrium solid {cof}, 
and how they change at the saddle point of the diffusional jump {uif} [127]. The 
expression for the temperature-independent prefactor of Eq. Ill was found to be 

e+ A5 V *B = , (112) 



where the 3N th frequency of the equilibrium crystal is associated with the attempt 
frequency of the diffusive jump, Tq. Schober, Petry and Trampenau assessed the 
vibrational entropy and enthalpy of vacancy migration using a lattice Green's func- 
tion method. For an idealized activation barrier of vacancy migration, they argued 
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that the different phonons should be weighted by the phonon DOS with a factor 
of l/(mw) 2 . This weighting of the phonon DOS accounts for how thermal energy 
causes larger mean-squared displacements for phonons of lower frequencies [128]. 
Defect formation entropies of order 2 k B /atom were reported, but these values are 
not quantitative. 

10.2.2. Solubility 

The thermodynamics of solubility can be understood with the help of Eq. 48, using 
Eqs. 49, 50, and 52. In the dilute limit of B-atoms in an A-matrix (small c), the free 
energy is 

F{c) ~ ziV/2[eAA + 2c(eAB -e AA )] 

+ k B TN[clnc + (1 - c)ln(l - c)] 

+ k B TVN[\nuAA + 2c\n(u AB /oJAA)] ■ (113) 

To calculate the solubility, we find the difference in free energy of a dilute A-B 
alloy with respect to pure A and pure B. 30 For illustration, assume all materials 
have the same crystal structure, and we know that the pure elements have zero 
configurational entropy of mixing. We calculate F(c=0) and F(c=l), weight them 
by their proportions, and take the difference 



AF(c) = F(c) - (1 - c)F(O) - cF(l) , (114) 

V^AA^BB / 

+ N k B T [c\n c + (1 - c) - c)] . (115) 



AF(c) = —czN2V + cDNk B T 



Seeking the minimum of AF(c) as dAF(c)/dc = 0, the equilibrium concentration 
for small c is 



2 i rr 

The exponential in Eq. 116 includes the interchange energy of Eq. 42, which is 
negative for unmixing systems. The prefactor is exp(— AS^/fee), comparable to 
Eq. 112. The configurational entropy does not appear explicitly in Eq. 116, but the 
derivation assumed that all sites of the A-crystal are equivalent for placing B-atom 
solutes. 

An analogous expression for a dilute concentration of B-atoms in an A-atom 
matrix, in equilibrium with a second chemically-ordered phase of precise composition 

Ai— jBj;, is 

/ F(B) - ±F(A^ x B x y 
c = exp 



exp 



k B T J ' 

Svib(B) — —S v it,(Ai- x B x ) 



exp 



E(B) - l£(Ai_ 
k B T 




30 The Hume-Rotherey rules, classic arguments for the size of {e}, and perhaps extendable to 
{lo}, are presented in Section 21.2. 
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where E(B) denotes the energy of formation for placing a B-atom on a matrix site 
and E(Ai- x TS x ) denotes the energy of formation per atom of the phase Ai-^B^ 
(similarly for F and SVib)- Any additional sources of entropy, such as electronic 
entropy, are added to SVibi except for configurational entropy, which was implicitly 
included when deriving Eqs. 116 and 117. 

Equation 118 can be used to obtain vibrational entropies if the energies of forma- 
tion are known. Because SVib vanishes at low temperatures, experimental approaches 
can isolate the E in Eq. 118. More typically, a plot of lnc versus l/T can be 
extrapolated to l/T = 0, giving the entropy difference of Eq. 118. Clarence Zener 
took this approach many years ago, and identified effects of vibrational entropy on 
solubilities in dilute Al alloys [129]. 

Section 1 discussed a key difference between vibrational entropy and configura- 
tional entropy that is usually used to calculate solubilities of dilute solutions. The 
vibrational entropy has an intrinsic temperature dependence, whereas the configu- 
rational entropy does not. The relative importance of vibrational entropy therefore 
increases with temperature, so it is sometimes important to account for a tempera- 
ture dependence in the difference «Svib(B) — ^S v ib(Ai- x B x ) of Eq. 118. For harmonic 
solids, this difference would be a constant at high temperatures. Nevertheless, an- 
harmonic effects can be large. It is plausible that changes with temperature of the 
vibrational entropy of mixing could cause a retrograde solubility. 31 

Huntington, Shirn and Wajda studied the vibrational entropy of formation of 
lattice defects [130]. They used an Einstein model with a Born-Mayer potential, 
and considered the changes in vibrational frequencies of individual atoms around 
a defect, or the frequencies of an interstital atom. These entropies were of order 
+1 fee/atom, with the defects giving the larger entropy. Freedman and Nowick [131] 
built on this work by analyzing data on solubility versus temperature for a set of 
metallic solid solutions to extract the "excess entropy," (i.e., entropy that was not 
accounted for by the entropy of mixing, Eq. 11). These contributions were generally 
positive, meaning that they added to the entropy of Eq. 11, but ranged from - 
1.1 &B/atom for Ni in Ag to +3.0 /ce/atom for Fe in Cu. Their analysis of vibrational 
entropy used the Einstein model, but it made the "rather drastic assumption" that 
there was no overall change in the vibrational frequencies of the solute atoms in the 
matrix, compared to the averages of the pure elements. This assumption permitted 
an analysis like that of Zener, which was based on the temperature dependence of 
the shear modulus. The effects of vibrational entropy seemed to be large. 

Desorbo [132] and Oriani [133] used calorimetric methods to estimate the excess 
entropy of near-equiatomic Au-Ni, and attributed it to vibrations. Desorbo found 
the excess entropy to be approximately +0.3 fce/atom at room temperature, and 
Oriani states, "The positive value of this deviation and the large magnitude of 
the vibrational entropy of solution are related to the loosening of the lattice that 
accompanies formation of the solid solution." 

31 Of course there is a drastic change in the temperature dependence of this difference at low 
temperatures below the Debye temperature (cf., Fig. 8b), but such low temperatures are often 
below the range of validity of T-c phase diagrams. 
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10.3. Carbides, Nitrides, and Oxides 

Fernandez Guillermet and Grimvall examined thermodynamic data on transition 
metal carbides, and reported some trends [134]. The standard entropy, primarily 
vibrational in these ordered structures, and the standard enthalpy were shown 
to depend on the electron to atom ratio after the carbides were separated into 
two groups - one with Sc, Ti, V, and the second with Cr, Mn, Fe, Co, Ni. The 
first group showed stronger variations in the standard enthalpy and entropy with 
electron-to-atom ratio than did the second group. For carbides in the first group, 
standard entropies of formation varied by typically 0.4 fce/atom for different carbon 
stoichiometrics, although the crystal structure also changed. Fernandez Guillermet 
and Du assessed the thermodynamic properties of the Fe-N system, and noted a 
large increase in Debye temperature as the nitrogen concentration was increased 
from the Fe4N to FeN3 (from 300 to 584 K) [135]. In the high temperature limit, 
this corresponds to a large difference in vibrational entropy of 3/cb ln(300/584) /atom 
= — 2.0 k~B /atom for these two compounds. 

The pressure-induced phase transformation of ZnO from the wurtzite structure 
to the rocksalt structure was the subject of a first-principles computational study 
using the planewave pseudopotential code VASP, with the code PHONON to cal- 
culate phonon frequencies by the direct method [136]. Vibrational entropy was 
calculated with the harmonic expressions for the phonon frequencies obtained at 
different volumes (i.e., the quasiharmonic approximation). The transition pressure 
was found to decrease with temperature (consistent with experimental results), and 
this was attributed primarily to a larger vibrational entropy of the high-pressure 
rocksalt structure. The ZnO-MgO composition-temperature phase diagram was also 
calculated, and it was found that vibrational entropy increases the solubility of ZnO 
in the rocksalt MgO by about a factor of two at a given temperature, whereas 
vibrational entropy decreases the solubility of MgO in wurtzite ZnO by a similar 
factor [136]. 

Some of the expectations and rules for metal systems may not be transferrable to 
nitrides or oxides. Structural units in oxide crystals, i.e., metal-oxygen octahedra, 
can undergo dynamics at low energies by rotating as a rigid unit, whereas the modes 
involving relative motions of the metal and oxygen atoms typically occur at much 
higher energies. These "rigid unit modes" of the octahedra [137] can break the 
symmetry of a structure of packed octahedra, and can be responsible for displacive 
phase transformations in minerals, for example [138]. In structures where some 
oxygen atoms at corners of the octahedra are not shared by other octahedra, there 
is more freedom in these reorientations. The material ZrW^Os has such a local 
structure, and it exhibits an unusual negative thermal expansion at low tempera- 
tures. An analogy is opening a pair of scissors, where the total length is reduced as 
the device is opened, but there is enough open space at the side to accommodate 
the increase in width. The relevant low energy vibrational modes in ZrW^Os have 
been identified by inelastic neutron scattering [139, 140], and they show a stiffening 
with increased temperature. A similar negative thermal expansion was observed for 
the related compound H1M02O8, and a similar interpretation was offered in terms 
of low-frequency optical modes [141]. For quite some time it has been known that 
Rbl exhibits a negative thermal expansion, as does silicon at low temperatures. For 
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Rbl the mechanical analogy is to a guitar string, where the increased vibrational 
amplitude of the string pulls inwards its points of suspension [139, 142]. 

11. First-Principles Studies of Solutions and Studies of Ordered 

Compounds 

11.1. Solubility in Al Alloys 

The interatomic force constants and the phonon DOS of fee Al metal were obtained 
by analyzing the thermal diffuse x-ray scattering from a single crystal of Al [143]. 32 
Neutron inelastic scattering experiments were later performed on Al, giving similar 
results for the lattice dynamics [144, 145]. Many features in the phonon dispersion 
curves were identified as Kohn anomalies from electron-phonon interactions (Sect. 
20). These anomalies are small and difficult to observe, but in principle they could 
have some thermodynamic effect. 

Ab-initio electronic structure calculations have proved valuable for calculating 
the energy and vibrational entropy of a matrix with solute atoms, and Eq. 118 
can be used to find the effects of vibrational entropy on solubility. This approach 
was used by Ozolins and Asta, who found a large effect of vibrational entropy on 
the solubility of Sc in Al [146, 147]. Although the solubility of Sc in Al is low, 
calculations without vibrational entropy predict it to be about 27 times lower than 
has been measured, but the inclusion of vibrational entropy brings good agreement 
with experient. The Sc dissolved in fee Al is in equilibrium with Sc in the LI2- 
ordered compound AI3SC, so the energies and entropies of ordered compounds were 
also calculated, as discussed below. The LI2 structure of AI3SC proved to be very 
stiff, with a high bulk modulus and small volume per Sc atom The phonon DOS 
of the LI2 structure had an average frequency higher than either fee Al or hep Sc, 
giving it a negative entropy of formation. This entropy of formation of AkjSc from fee 
Al and hep Sc was big, — 0.70 fee/atom. The strong hybridization between Sc d- and 
Al p-electrons is responsible for short and stiff bonds in this intermetallic compound. 
For Sc solutes in fee Al, stiffer first-nearest-neighbor bonds are also expected, and do 
occur. Nevertheless, the vibrational entropy per Sc atom in fee Al is approximately 
+0.4 fee/atom, indicating an overall softening of phonon frequencies. The softening 
occurs for Al-Al bonds in the more distant neighbor shells around the Sc atom. Two 
inter-related contributions were identified as a source of this softening: 

— The element Sc is larger than the Al atom, and tends to increase the lat- 
tice parameter of the alloy. The elongation of the Al-Al bonds is a source of 
softening. 

— Al is more electronegative than Sc, so we expect a charge transfer and a charge 
perturbation around Sc atoms. The conduction electrons of Al redistribute 
around the charge perturbation, much as Friedel oscillations in a free electron 
gas. The entropic effects of the charge disturbances are of similar size as the 
effects of volume expansion. 



This was sophisticated work, especially in its day, and few other phonon studies have since 
been done this way. 
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Ozolins, Sadigh, and Asta [148] performed a computational study on the sol- 
ubility of Si in Al. The Al-Si phase diagram is a simple eutectic, so the relevant 
phases are fee Al, fee Al-Si, and dc Si. The phonons in dc Si are quite stiff in 
comparison to the other phases. Removing a Si atom from covalently-bonded dc Si 
and placing it in fee Al has the net effect of lowering the vibrational frequencies of 
all atoms. The vibrational frequencies at Al atoms in the vicinity of Si solutes are 
not altered significantly. (There is no large charge redistribution as for the case of 
Sc, for example.) With the dominant effect of solubility being the loss of short, stiff 
covalent bonds in dc Si, the vibrational entropy of formation for dilute Si in Al was 
calculated as being +2.6/cB/Si atom, quite a large effect. This seems qualitatively 
consistent with the softening of transverse acoustic modes measured by inelastic 
neutron scattering on supersaturated Al-Si solid solutions [149]. It should be noted 
that Ozolins, et al., encountered difficulties with the LDA approximation in this 
study of Si solubility, owing to differences in the exchange-correlation energy caused 
by charge inhomogeneities in the different phases. 

11.2. Aluminum Compounds 

Aluminum alloys are of high technological importance because pure Al is mechani- 
cally soft, and second phases can provide strengthening mechanisms. One of the best 
known examples of precipitation hardening occurs for Al-Cu with the nonequilibrium 
phases 9" and 9'. These eventually transform to the equilibrium 9 phase, so the 
free energy of the 9 phase must be lower than that of the metastable 9' phase. 
Nevertheless, Wolverton and Ozolins calculated the electronic energy of the AI2CU 
9 phase, and found it to be 15meV/atom larger than for the 9' phase [150]. The 
stability of the 9 phase was attributed to its vibrational entropy, approximately 
+0.37 &b /atom larger than that of 9' in the harmonic approximation. An entropic 
stabilization of the 9 phase implies a phase transformation to the lower-energy 
9' phase upon cooling. The transformation temperature is predicted to occur at 
T c = 15 meV/0.37 k% ~ 200°C. There is weak experimental evidence for this phase 
transformation, but it should be quite sluggish at this low temperature. 

The entropic stabilization of the 9 phase may originate with features of its low- 
symmetry structure. The lowest-frequency optical modes in the 9 phase allow for 
perpendicular motions of neighboring Al and Cu atoms, with minimal change in 
Al-Cu bond lengths. Wolverton and Ozolins suggest that vibrational entropy stabi- 
lization may be important for other low-symmetry structures such as Laves phases 
in transition metal systems. The general idea is that the low symmetry may result 
in a less favorable bonding energy, but this is compensated at higher temperatures 
by a larger vibrational entropy originating from floppy atom motions. 

Vibrational entropies of formation were calculated for several intermetallic com- 
pounds in the Al-Sc system by Asta and Ozolins [147]. Although the vibrational 
entropy of formation of the stiff LI2 AI3SC from fee Al and hep Sc had the largest 
magnitude, —0.77 fee /atom, the B2 AISc and orthorhombic AISc compounds had 
vibrational entropies of formation of — 0.58 &b /atom and — 0.24 k-Q /atom. The or- 
thorhombic Bj compound (BCr prototype) has a lower symmetry than the cubic 
LI2 and B2 structures, and its low-frequency optical modes give it the least negative 
vibrational entropy of formation of these compounds with strong Al-Sc bonds. 
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Figure 12. Vibrational entropies of mixing of Hf in Al. The symbols are for ordered structures. 
Open squares were calculated by the transferable force constant approach while the solid circles 
are the values predicted by the cluster expansion method The solid lines are predictions from the 
cluster expansion method for a disordered fee solid solution of Al-Hf. After [45] . 

Quasicrystals, or icosahedral-phase compounds, are found in a number of aluminum- 
transition metal systems as metastable phases, and are in fact stable phases in some 
systems such as Ai-Co-Ni. In the case of AlCo^, where 0.23 < x < 0.26, a large 
number of crystal structures have been reported, some of which are quasicrystal ap- 
proximants. It is challenging to determine which of these is the thermodynamically- 
stable structure. Energetic comparisons with ab-initio calculations can help. In such 
work by Mihalkovic and Widom [151], the orthorhombic and monoclinic AI13C04 
structures were found to be unstable energetically. Nevertheless, a combination of 
defect entropy and vibrational entropy made these AI13C04 structures stable at 
temperatures above 750 K. 

As discussed in Sect. 5.3, Liu, et al., identified chemical modifications that im- 
prove the agreement between the bond-stiffness-versus-bond-length model and first 
principles calculations for Al-TM alloys, where TM = {Ti, Zr, Hf} [45]. It is possible 
that these concentration dependences are correlated to the difference in electroneg- 
ativity between the TM and the Al atom - this difference is smallest for Al-Ti, 
as are the concentration dependences of the curves of bond stiffness versus bond 
length. Consider, for example, the systematics of Al-Hf shown in Fig. 12. The strong 
hybridization of Al p- and Hf d-electrons tends to build up electron density in the 
bonding levels between Al-Hf INN pairs, causing a stiffening of phonons and a 
negative entropy of mixing. Figure 12 shows that the entropies of mixing are large 
and negative in all cases, typically -1 fee /atom. The systematics of ordering are 
less obvious. The ordered phases tend to be stiffer, with negative vibrational en- 
tropies of ordering of approximately -0.2 fee /atom, but there are exceptions, perhaps 
associated with lower symmetry structures. 

11.3. Low-Symmetry Phases at High Temperatures 



Any survey of alloy phase diagrams shows examples of chemically-ordered com- 
pounds that exist at high temperatures, but decompose to other phases at lower 
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temperature. High-temperature phases are stabilized by entropy, but there is no 
configurational entropy for chemically-ordered structures, of course, and vibrational 
entropy likely plays a role in their high-temperature stability. 33 Many of these high- 
temperature ordered phases are low-symmetry crystal structures. The B81 NiAs 
structure is such an example, and may sometimes be stabilized by vibrational 
entropy. The B81 structure has a number of asymmetrical intersticies that could 
accommodate large-amplitude atom vibrations, which may be of low frequency 
and anharmonic. More interestingly, in about half the phase diagrams where the 
B81 structure appears, it is not stable at low temperatures, but only appears at 
high temperature. Vibrational entropy seems a likely source of stabilizing the B81 
structure in AlCu, CrSe, CuSn, FeS, Culn, InPt, PdSb, FeSn, CrTe, FeTe, MS, CrS, 
CoSe and FeSe, which appear only at high temperatures [152, 153]. Many of these 
compounds exist off the equiatomic stoichiometry, however, suggesting there is also 
a role for defect entropy in their high-temperature stability. 

At high temperatures, the full thermodynamics is likely more complex. For exam- 
ple, the B8i NiAs structure of CuSn undergoes an ordering transformation to a long- 
period structure at temperatures below 187° C, with an ordering enthalpy measured 
to be —0.26 ± 0.03 kJ/mole [154], implying an ordering entropy of — 0.07 &b /atom. 
This is an unusually small entropy for chemical ordering, suggesting that the low- 
temperature ordered phase has a relatively large vibrational entropy (or only a small 
fraction of atoms are involved in the ordering process, implying that the chemical 
disorder at high temperatures is small). The B8i phase of FeSe is a case where a 
number of other sources of entropy may affect its stability. 34 It is an antiferromagnet 
with a Neel temperature of 600 K, which is very close to the temperature where the 
B8i phase appears on the phase diagram (633 K). Furthermore, the Fermi level of 
B8i FeSe lies in a deep valley between two closely-spaced, high peaks in the electronic 
DOS [157]. Thermal smearing, discussed below in Section 20, should increase the 
DOS at the Fermi level, enhancing the electron-phonon coupling, although this trend 
should reverse at higher temperatures. Both magnetism and the adiabatic electron- 
phonon interaction should be considered when assessing the stability of the B8i 
phase of FeSe. 

As a curious counterexample, a high-symmetry B2-ordered phase exists in the 
Cu-Zr system only above 715°C at 51.5% Zr [152]. This B2 phase undergoes a 
eutectoid transformation below 715°C to two other intermetallic compounds of lower 
symmetry [158, 159]. Thermodynamic assessments give it an entropy of stabilization 
of 0.46 fee /atom with respect to pure Cu and Zr [160], which seems too large to be 
accounted for by configurational entropy (unless there is a high degree of chemical 
disorder in the B2 structure). Although it seems reasonable to propose a rule of 
thumb that vibrational entropy helps stabilize low-symmetry structures at high 
temperatures, the B2 phase in Cu-Zr serves as a warning. Another example is the 
hexagonal C32 structure. It is not a low-symmetry structure, but many C32 diboride 

33 Some partial chemical disorder may contribute configurational entropy, and antisite defects 
may also influence the phonon dynamics in ways that are not understood. 

34 It has recently been found that the B10 structure (PbO prototype) of Fe-Se, stable at more 
Fe-rich compositions, is a superconductor with similarities to the new iron arsenates [155, 156]. 
Perhaps other charge or spin excitations in the B8i phase of Fe-Se are relevant for phase stability. 
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phases are stable to extremely high temperatures. Some, like M11B2, are unstable at 
lower temperatures. 

12. Chemical Order-Disorder Transformations 

12.1. Transition Metal Aluminides 
12.1.1. Fe 3 Al 

When the atom displacements, u, are comparable to the inverse of the scattering 
wavevector, Q, an elastic Bragg diffraction is suppressed by the Debye- Waller fac- 
tor exp (— Q 2 (u 2 }), as discussed in Sect. 8.1. A closely-related effect occurs with 
the local probe methods "extended x-ray absorption fine structure" (EXAFS) or 
"extended electron energy- loss fine structure" (EXELFS), which measure scattering 
above a core ionization edge of an atom. The electron leaving the ionized atom has a 
wavelength that allows self-interference by backscattering from the nearest-neighbor 
shells of atoms. For energies E above the absorption edge E a , the wavevector, k, of 
the outgoing electron is E — E a = h 2 k 2 (2m e )~ 1 . When 2n/k matches the roundtrip 
nearest-neighbor distance, constructive interference occurs, and the ionization prob- 
ability is slightly enhanced. 35 The Fourier transform of the EXAFS or EXELFS 
oscillations above E a can be used to identify nearest-neighbor distances and the 
scattering strengths of the neighboring atoms. If the nearest-neighbor atoms are 
vibrating with respect to the ionized atom, the self-interference of the outgoing 
electron is suppressed. The effect is essentially the same as a Debye-Waller factor, 
but it is associated with the mean-squared relative displacement (MSRD) between 
neighboring atoms. Measurements of the MSRD as functions of temperature for Fe 
and Al atoms in DC>3-ordered and disordered solid solutions of FesAl are shown 
in Fig. 13 from work in 1990 [161]. What is important is that the MSRD of Al 
atoms is larger, and more temperature-dependent, in the disordered than in the 
ordered phase. A similar but weaker trend is found for Fe atoms. Fitting these data 
to MSRD curves calculated with a local Einstein model provides characteristic tem- 
peratures listed in the caption. With these Einstein temperatures we can estimate 
a difference in vibrational entropy of 0.4 ke/atom, with the disordered material 
having the larger entropy. These EXELFS measurements were repeated with better 
experimental technique, but with similar results [162]. 

Subsequent calorimetry work using Eq. 110 [163] showed that the difference in 
vibrational entropy of disordered and ordered FesAl is (0.1±0.03) ke/atom, much 
smaller than the 0.4 ke/atom deduced from the simple EXELFS interpretation 
with an Einstein model. The discrepancy gives insight into the origin of the change 
in vibrational entropy. The EXELFS method is more sensitive to high-frequency 
longitudinal phonons than other phonons [164]. If disordering affects primarily the 
high-frequency phonons, EXELFS will tend to overestimate the change in vibrational 
entropy. 

Figure 14 presents results on 57 FesAl from nuclear resonant inelastic x-ray scatter- 
ing (NRIXS) measurements, and partial DOS curves calculated for the D03-ordered 



There are numerous other factors to consider, especially the change in phase of the electron 
wavefunction upon scattering, but these are beyond the scope of the present article. 
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Figure 13. (a) Experimental MSRD between Al atoms and their Inn shell relative to data at 96 K. 
Lines are calculations of MSRD based on Einstein temperatures of 350 K for Al in disordered Fe3 Al, 
and 460 K for Al in D03-ordered FesAl. (b) Experimental MSRD between Fe atoms and their Inn 
shell relative to data at 96 K. Lines are calculations of MSRD based on Einstein temperatures of 
400 K for Fe in disordered Fe 3 Al, and 440 K for Fe in D0 3 -ordered Fe 3 Al. 

structure [165]. The NRIXS method is similar to incoherent inelastic neutron scat- 
tering, but the spectral intensities of NRIXS are weighted entirely on the 57 Fe atoms, 
which exhibit the nuclear resonance. The calculated partial phonon DOS of Al shows 
that the optical modes at energies around 43 meV are dominated by the motions of 
the Al atoms. The shapes of the Fe and Al partial DOS curves at the energies of the 
optical modes are rather similar. 36 The 8[c] Fe atoms, which are first neighbors of the 
Al atoms, have small displacements in response to the Al vibrations. Nevertheless, 
the partial DOS of 57 Fe is measured reliably owing to the low background and 
high sensitivity of the NRIXS method. Figure 14 shows that with disordering, 
the spectral contribution from the optical modes becomes less distinct, but more 
importantly these modes shift to lower energy. This effect depends primarily on the 
short-range order in the alloy, and substantial changes to the phonons occur well 
before any superlattice diffractions from long-range order are observed in diffraction 
measurements [166]. A downwards shift of about 10% for a quarter of the modes 
accounts for a change in vibrational entropy of —0.25 x 31n(0.9) = 0.08 fee /atom, in 
good agreement with the calorimetry result of (0.1±0.03) ke/atom. It is unfortunate 
that the partial phonon DOS of Al cannot be measured, but the MSRD data of Fig. 
13a also point to changes in the Al atom dynamics at high frequency as a main 
source of the vibrational entropy of the order-disorder transformation in FesAl. 

12.1.2. Ni 3 Al 

Section 4.2 discussed how phonon frequencies are expected to increase with the 
strengths of chemical bonds, or pairwise interaction energies. For an alloy with 
strong ordering tendencies (a large V and a high critical temperature (Eqs. 62, 64, 
70)), we might expect a large change in vibrational entropy upon disordering. The 
critical temperature for DO3 ordering in FesAl is 820 K, although B2 order persists 
to 1100 K or so. On the other hand, for M3AI the LI2 ordered structure persists to 
melting at 1635 K, and chemical disorder is notoriously difficult to induce in N13AI. 

36 Three phonon branches out of twelve (4-atom unit cell) correspond to these optical modes, so 
from 40-45 meV the partial DOS of Al must sum with the Fe partial DOS to account for 25% of 
the total DOS. 
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Figure 14- Phonon partial DOS curves for 57 Fe in 57 Fe3Al. Top: Calculated with a Born - von 
Karman model using force constants from prior neutron scattering measurements on single crystals 
of DC>3-ordered FeaAl, convoluted with an experimental resolution function. Bottom: 57 Fe phonon 
partial DOS curves from FeaAl with high long-range order ("Ord"), some short-range order ("Par"), 
and no measurable short-range order ("Dis"). 



Disordered samples of Ni-25at.%Al were prepared by physical vapor deposition onto 
cold substrates, and their heat capacity was compared to annealed samples with LI2 
chemical order using differential scanning calorimetry (Sect. 7) and Eq. 110 [167, 
168]. X-ray Debye- Waller factor measurements and EXELFS MSRD measurements 
were also performed on these pairs of ordered and disordered samples, indicating 
that the vibrational entropy of the disordered Ni 3 Al was 0.3 fee/atom larger than 
the ordered (error bars are probably ±0.1). 

Inelastic neutron scattering measurements were performed on disordered N13AI 
prepared by high-energy ball milling, and experimental spectra are presented in Fig. 
15a [169]. The optical modes at 39meV in the ordered alloy are clearly suppressed 
by disordering, and the disordered spectra show new intensity around 32 meV. The 
phonon DOS of the ordered N13AI was possible to deduce from the experimental 
spectra by a neutron weight correction based on a comparison with a Born-von 
Karman calculation from prior work on single crystals of LI2 M3AI [170]. As ex- 
pected by footnote 36 for ordered F83AI, there is indeed about 25% of the intensity 
in the optical mode peak around 39 meV in the neutron-weight corrected spectrum 
in Fig. 15b. The proper correction for the disordered sample is less obvious, but it 
should lie between two extremes. The first, labeled "Dis 1" in Fig. 15b, assumes a 
virtual crystal, with no neutron weighting of the different parts of the phonon DOS. 
The intensity at high energies is probably too low, since many Al atoms are still 
expected to vibrate in local modes. The second, labeled "Dis 2" in Fig. 15b, used the 
same neutron weight correction as for the ordered alloy. This method probably gives 
too much intensity in the modes around 39 meV because with disorder there is clearly 
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Figure 15. (a) Neutron energy loss spectra for ordered and disordered Ni 3 Al powders, obtained for 
four values of Q. Data for the ordered alloy are offset by 1500 counts, (b) Phonon DOS curves of 
ordered (LI2) and disordered N13AI (Disl, Dis 2). Neutron weight corrections were performed for 
the ordered data. There is uncertainty in the neutron weight correction for the disordered alloy, 
with two extremes shown as Dis 1 and Dis 2. 



some intensity that shifted to lower energies, which are suppressed somewhat by the 
renormalization of the DOS curve to unit area. The effects of disorder on the optical 
modes of N13AI, which are dominated by the motions of Al atoms, are larger than 
the similar effects found for FesAl, consistent with the calorimetric measurements 
of a vibrational entropy of ordering that are a factor of 2 or 3 larger for M3AI than 
Fe 3 Al. 

After the first report that vibrational entropy was important to the thermody- 
namics of ordering in N13AI [167], a number of computational studies addressed the 
entropy difference of ordered and disordered Ni 3 Al [171, 172, 173, 174, 175, 176]. 
Agreement with experiment, or with each other, was not good. A density functional 
theory study on M3AI by van de Walle, Ceder and Waghmare [176], arguably the 
most sophisticated at its time, found essentially no difference of vibrational entropy 
between the ordered and disordered materials. The biggest effects were in the optical 
modes involving Al atom motions, as was found experimentally, but in the DFT cal- 
culations these modes were shifted to higher energies in the disordered material. One 
source of this discrepancy is the difference in lattice parameter of the ordered and 
disordered N13AI. Phonon frequencies are highly sensitive to lattice parameter, as 
discussed in Sect. 5.2. Experimentally, the lattice parameter of the ordered material 
was larger than the disordered, but the DFT calculations gave the opposite result, 
and EAM calculations gave an even larger difference [172], again in the opposite di- 
rection as experiment. With Eq. 68 we obtain 5S = 31n(l+7(5V /Vq) = 0.17 k^/atom., 
using 5V/Vq = 0.03 for the difference between ordered and disordered material pre- 
pared by high-energy ball milling, and assuming 7 = 2. Unfortunately, this worsens 
the agreement between theory and experiment, so a more plausible reason is that the 
chemical disorder in ball milled M3AI is not of the antisite variety assumed in the 
computational work. This makes it difficult to do a rigorous comparison between 
theory and experiment. Another peculiarity of the experimental phonon DOS of 
the disordered material in Fig. 15b is the enhanced intensity at energies below 15 
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meV. This is an interesting aspect of nanostructured materials, discussed below in 
Section 16.2. The positions of the Van Hove singularities in the phonon DOS of 
nanostructured materials are usually the same as for bulk materials, however, so 
the effects on optical modes, other than broadening, are not expected to be large. 

12.2. B2 Ordering in CuZn and FeCo 

In 1938, a few years after Bragg and Williams [177] and Bethe [178] explained order- 
disorder transformations without vibrational entropy, Nix and Shockley argued that 
the vibrational entropy change of order-disorder transformations could be assessed 
from Debye temperatures [179]. Zener assessed the temperature-dependence of the 
phase boundaries of disordered a- and /3-phase Cu-Zn, and argued that the (3 
phase must be stabilized in part by vibrational entropy (which works out to be 
approximately 0.13 k^/atom) [180]. The role of vibrational entropy in the chemical 
order-disorder transformation of /3-brass was considered with increasing sophistica- 
tion in the mid-twentieth century. In 1953 Sarolea used a Lennard-Jones model and 
a method of physical chemistry to deduce that vibrational entropy would cause a 
small 4.5% decrease in the ordering transition temperature in CuZn [181]. In 1955, 
Booth and Rowlinson developed a classical partition function for vibrations in the 
Einstein model to extend the Bragg- Williams model of ordering. They estimated a 
50% increase in the heat capacity of /3-brass across the order-disorder transition, 
and a lowering of the critical temperature [182]. Using a specialized calculation, 
Wojtowicz and Kirkwood found that the decrease in critical temperature was about 
a factor of two [183]. The B2 ordering transformation in Cu-Zn seems worthy of 
study by more modern methods. 

The alloy Fe-50at.%Co has a B2 order-disorder transformation at 1000 K. The en- 
thalpy of the ordering transformation, Ai? or( j, was measured by differential scanning 
calorimetry, and found to be 4.4kJ/mol [184, 185] giving a transformation entropy 
A>Sord = Ai/ or( j/T or d = 0.53 fee/atom. This could be accounted for by the configura- 
tional entropy alone, but other sources should be assessed. Both the disordered and 
ordered phases are ferromagnetic at T orc j. It seems that the magnetic entropy does 
not change significantly with ordering, since the saturation magnetization has little 
change around 1000 K [186], and in-situ Mossbauer spectrometry studies show that 
the ordering process at 800 K does not cause any discontinuity in the temperature 
dependence of the hyperfine magnetic field [187]. Although no phonon measurements 
have been performed at elevated temperatures, phonon DOS curves of ordered and 
disordered Fe-50at.%Co alloys have been measured at 300 K [188]. The features of 
the phonon DOS curves are somewhat sharper in the B2-ordered sample, but there is 
no net shift of the phonon frequencies that could account for a significant difference 
in vibrational entropy between the ordered and disordered materials. Interestingly, 
this is qualitatively consistent with the prediction of Eq. 70 - unlike the general 
trend for ordered compounds, the B2 ordered phase of FeCo has a lattice parameter 
that is 0.2% larger than for the bcc solid solution [189]. With Eq. 69 this dilation 
predicts a small entropic contribution of perhaps 0.03-0.04 fee /atom. The role of 
vibrational entropy on the B2 ordering transformation in FeCo appears to be small. 
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12.3. Ll 2 Ordering 
12.3.1. Cu 3 Au 

In 1954, Bowen reported measurements of the electrical resistivity of ordered and 
disordered CU3AU over a range of cryogenic temperatures [190]. Debye temperatures 
of 175 and 197 K were reported for the disordered and ordered alloys, suggest- 
ing that the disordered phase had a larger vibrational entropy by the amount 37 
3 ln(197/175) = 0.35 [k B /atom]. Flinn, McManus and Rayne [191] measured the 
elastic constants of ordered and disordered CU3AU from 4.2 to 300 K, finding dif- 
ferences of less than 1%. They criticized Bowen's results, "It is concluded that 
little significance can be attached to [Debye temperatures] obtained from resistance 
measurements on CU3AU." 38 Rayne had performed an earlier calorimetric study on 
ordered and disordered CU3AU [193], and found evidence that vibrational entropy 
of ordering was between the values reported in [190] and [191]. Nevertheless, the 
message from Flinn, McManus and Rayne was that the difference in vibrational 
entropy of ordered and disordered CU3AU was not important thermodynamically. 
This thinking persisted for decades. 

Differences in the phonon spectra of ordered and disordered CU3AU were mea- 
sured more recently [194] with a triple-axis spectrometer, but the thermodynamic 
consequences were not assessed. There is risk, however, in interpreting the peaks in 
the phonon scans at constant Q, measured on single crystals of disordered material. 
Specifically, the peaks in 1(e) from the disordered solid solution of CU3AU were inter- 
preted as energies of phonons in a "virtual" fee crystal, which assumes one species 
of atom with average properties of Cu and Au atoms. The virtual crystal is not 
expected to be trustworthy when the Cu and Au atoms differ in mass by a factor of 
three, and differ considerably in their atomic radius. Nevertheless, this approach was 
used in an earlier reconciliation of the phonon dispersion curves with heat capacity 
data [168, 195]. The calorimetric measurements by Nagel, et al., showed that as a 
disordered solid solution, CU3AU has a larger vibrational entropy than for the LI2 
structure by about 0.14±0.05 /ce/atom [195, 196]. A similar value was obtained from 
a lattice dynamics calculation with forces derived from tight-binding potentials with 
a virtual crystal approximation [197] and other potentials [174, 198], although there 
were differences in interpretation. Surface effects were also investigated [199]. 

The neutron weighting problem (Sect. 8.3) is a major challenge for phonon spectra 
of CU3AU measured by time-of-flight instruments. An atom of Cu is 3.2 times more 
efficient than Au for phonon scattering. Because the atom masses differ by a factor 
of 3.1, we expect a tendency of Au atoms to have larger displacements in low- 
energy modes, and Cu atom motions to dominate at higher energies. Correcting for 
neutron weighting is not difficult when a Born-von Karman model is developed for 
ordered phases with a periodic unit cell, but it is difficult to perform for disordered 
solid solutions. It is probably more reliable to work with the entropy as an integral 
quantity of the phonon DOS, defining a neutron-weighted vibrational entropy, S^, 

37 The factor of 3 arises because each atom has 3 degrees of freedom in 3 dimensions. 

38 Rossiter and Bykovec later showed that the electrical resitivity of CU3AU depends in a complex 
way on the presence of antiphase domain boundaries, and is not a simple function of order parameter 
[192]. 



58 



B. Pultz 




that differs from the true vibrational entropy S. 



vib- 



cnw _ 3.2(1 x)S^(x) + xS^(x) moA 
6vib - 3.2(1 -x)+x ' (119) 

S vih = (1 - x)^(x) + xS^(x) , (120) 

where S^(x) and S^{J are the phonon partial DOS functions for Cu and Au at 
the Au concentration, x. 39 Several approaches were used to correct Eq. 119 for the 
neutron weighting of Cu-Au alloys [200], with the best estimate that the disordered 
phase has a larger vibrational entropy than the ordered by 0.15±0.05 fee/atom at 
300 K. 

A method for neutron weight correction of solid solutions is described in section 
13.2, and an early version of this approach was used to correct the phonon DOS 
curves in Fig. 16 [188, 200, 201]. The most prominent effect of disorder is a shift of 
intensity from peaks at 21 and 25 meV to lower energies, increasing the vibrational 
entropy of disordered CU3AU by about 0.15 k-Q /atom, in general agreement with 
the calorimetry results [195]. Other phonon spectra [200], and computational work 
by Ozolins, Wolverton and Zunger [202], also showed this change in the spectrum 
at 21meV owing to changes optical modes (or local oscillations) emphasizing Cu 
atoms. 

Using the force constants of Ll2-ordered CU3AU, the motions of atoms in the 
different modes can be identified. The highest frequency modes around 25meV 
involve the motions of Cu atoms against Au atoms. The Cu-atom motions were 
found primarily in the modes around 21 meV, and are associated with local Cu atom 
movements with less response from each Au neighbor. In some of these modes, a Cu 
atom moves perpendicular to the face of the of the standard LI2 structural cube. 
The constraints on the motions of Cu atoms from the framework of the massive 
Au atoms, a characteristic of the LI2 structure of CU3AU, could be removed in the 
disordered alloy. This may cause a shift of the modes at 21meV to lower energy. 
Although we do not have a definitive picture of the lattice dynamics of the disordered 



39 Note that when S°g(x) = S^(x), S vih = S% 
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alloys, the phonon DOS curves shown in Fig. 16 are consistent with this change in 
constraints of Cu atom motions. 

A resonance mode associated with the damped motions of the massive Au atoms 
is a prominent feature of the phonon DOS of disordered Cu-Au solid solutions with 
low Au concentrations. The broad resonance mode is centered at approximately 
9meV, and dominates the Au partial phonon DOS over a range of compositions 
[201, 203]. The resonance mode seems responsible for a more rapid rise in vibrational 
entropy as Au is added to Cu, as compared to adding Cu to Au. Unfortunately the 
study of Au-rich alloys by neutron scattering is difficult owing to absorption, and 
this work is incomplete. 

Finally, we return to the observations of Flinn, McManus, and Rayne [191], who 
measured no difference in the elastic constants of ordered and disordered CU3AU. 
From both theoretical and experimental lattice dynamics studies, it is also found 
that ordered and disordered CU3AU have minimal differences in the low-frequency 
phonons that determine the elastic constants. The higher frequency phonons are 
responsible for the larger vibrational entropy of disordered CU3AU. Clearly it is risky 
to use Debye temperatures from elastic constant measurements to assess differences 
in vibrational entropy of alloy phases. 

12.3.2. Pt 3 Fe 

To further explore the effects of local environments on the phonon DOS of the 
disorder-Ll2 order transition, a nuclear resonant inelastic x-ray scattering (NRIXS) 
experiment was performed on samples of Pt| 7 Fe [204] to obtain a phonon partial 
DOS of Fe atoms. In the LI2 structure of PtaFe, all 12 ffist-nearest-neighbors 
surrounding a relatively light Fe atom are massive Pt atoms. The structural re- 
arrangements of chemical disordering will vandalize this framework around the 
Fe atoms, so their vibrational frequencies are expected to change. In Ll2-ordered 
PtsFe, the phonon partial DOS of Fe atoms is largely a single peak at 25meV, 
approximately an Einstein mode. Figure 17 shows that with disorder, approximately 
10% of the intensity of this peak is transferred to a broader peak around 19 meV. 
This transfer of intensity depends primarily on the short-range order in the alloy, 
suggesting that local environments dominate the changes in the phonon partial DOS. 

A first interpretation of the spectra in Fig. 17 is that there are two vibrational 
environments for 57 Fe. There is one subspectrum at 25 meV associated with a nearly- 
full complement of Pt ffist-nearest-neighbors, and a second subspectrum at 19meV 
associated with enough Fe first-neighbors to vandalize the cage of Pt, and shift the 
vibrational frequency of the central Fe atom. If the disordered solid solution is a 
random solid solution, the numbers of different Fe environments can be calculated 
with a binomial probability distribution. From the areas of the two vibrational sub- 
spectra, the vandalized cage occurs with approximately 5 or more nearest neighbors. 
Of course the actual situation with the local vibrations is expected to involve the 
mutual positions of the different Pt atoms about the Fe atom, not just the total 
number of Pt atoms in a first-neighbor shell. Nevertheless, the results of Fig. 17 
suggest that the concept of a cage structure is relatively robust, and more than 
one or two misplaced atoms are required for a significant change in the vibrational 
frequencies of the Fe atoms. 
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Figure 17. Phonon partial DOS curves for 57 Fe in PtjpFe with LI2 order (ordered), with short-range 
order but little long-range order (partial), and prepared by piston-anvil quenching without any 
measurable short-range order (disordered) [204] . Also shown is a curve calculated from a Born-von 
Karman model (B-vK, which extends vertically to 1.8) [205]. 

12.3.3. Formation Entropy of LI2 Compounds 

Further information about the vibrational entropy of LI2 compounds is obtained 
by comparing the vibrational entropy of the compound to the vibrational entropies 
of the pure elements from which it is formed. The advantage of this comparison 
is that the phonon DOS curves of the pure elements are usually available, and 
are generally accurate. The systematics of the calculated "entropy of formation," 
AS*'?™, gives insight into the ordered structure, which helps our understanding of 
the order-disorder transformation itself. 

From phonon DOS curves published in the literature, or from DOS curves ob- 
tained from a Born-von Karman model with published force constants, the high- 
temperature limit of the vibrational entropy was obtained with Eq. 16 [206]. In this 
expression, the DOS g a (e) was the DOS of the ordered LI2 compound, with com- 
position A0.75B0.25, and gP(e) was the weighted sum of the DOS curves of elements 
A and B, 0.75g A (e) + 0.25<7 B (e). A number of correlations were explored between 
AS'^'™ and the properties of the constitutent atoms [206], but the correlation shown 
in Fig. 18 was most successful. 

The lattice parameter of the standard cube of the LI2 structure is usually similar 
to the lattice parameter of the fee unit cell of the majority species. Figure 18 
shows that when the minority species is larger than the majority species, there 
is a significant decrease in the vibrational entropy of formation. This originates with 
an overall stiffening of phonons, i.e., a shift of the LI2 spectrum to higher energies 
with respect to the average of the pure elements. On the other hand, inserting a 
smaller minority species between large majority atoms seems to cause only a modest 
softening of phonons in the LI2 structure, and this is not such a strong function of 
the difference in metallic radius. The systematics of Fig. 18 support the concept of a 
caged minority species in the LI2 structure - when the minority species is large and 
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Figure 18. Vibrational entropies of formation of LI2 alloys as a function of difference in metallic 
radius between the majority species and the minority species. The difference Ar mcta i is positive 
when the majority species is bigger than the minority species [206]. Line is from the second term 
of Eq. 124. 



"squeezed," the phonon frequencies are increased. An intuitive picture is that the 
atoms are slightly-compressible spheres, which have higher interatomic forces when 
they are compressed substantially 

The extended bond proportion model of Section 5.2 gives a semiquantitative 
interpretation of the trend in Fig. 18. In the LI2 structure, half of the INN pairs are 
A-B, half are A- A, and there are no B-B INN pairs. Because the lattice parameters 
of the LI2 compounds are similar to those of fee A, we assume that the INN A-A 
bonds between majority species are unaffected by the minority solute atoms (and 
we need not consider changes in the B-B bond length, only elimination of B-B bonds 
during LI2 compound formation). The vibrational entropy of formation per atom, 
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where the superscript "0" denotes the value for the pure elements or for Ar meta i = 0. 
The AV/Vo of Eq. 123 can be obtained approximately from the Ar meta i of Fig. 18 
as AV/Vo — 3Ar meta i/2ro, where 2ro ~ 3A, and assuming 7 = 2 
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The second term in Eq. 124 gives the dependence on Ar mcta i, and is shown as the 
diagonal line in Fig. 18. The first term in Eq. 124, the prediction from the bond 
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proportion model (cf., last term in Eq. 61), should cause the line to be offset from 
the origin. This first term is also expected to cause deviations from the line. 

12.3.4. Mass and Volume Effects in Intermetallic Compounds 

It is interesting to compare changes in the phonon DOS when one of the atomic 
species in an ordered compound is replaced with atoms of a different element. 
After this substitution, the vibrational entropy is expected to decrease with force 
constant, increase with mass, and increase with lattice parameter. All such effects 
could be present simultaneously, and in different amounts for different alloy systems, 
so rigorous predictions are challenging. Experimental phonon spectra were used 
to test these expected correlations when the minority atoms were replaced and 
the crystal structure was largely unchanged [206]. Changes in vibrational entropy 
were obtained when the minority B-atoms in A3B LI2, A3B A15, and A2B Laves 
phase compounds were replaced by a different species. The ASVib showed a fair 
correlation to the expected effects of a change in the mass of the minority species 
(3/2 c fee hi(mi/m2)), when added to the effects expected from the change in unit 
cell volume (cf., Eq. 66). A precise correlation is not expected, however, in part 
because the interatomic distances do not track the change in volume of the unit cell. 

An example of such a change in vibrational entropy with minority species sub- 
stitution is shown in Fig. 19 [207]. Figures 19a and 19b show the 57 Fe phonon 
partial DOS of Ll2-ordered Pt| Fe and Pd| 7 Fe, together with curves calculated 
with force constants from prior inelastic neutron scattering experiments [205, 208]. 
The atoms Pt and Pd have many chemical similarities, so it might be expected that 
the interatomic forces do not change significantly with the substitution of Pd for 
Pt. Figure 19c shows the result of a Born-von Karman calculation that used the 
force constants for Pt, but the mass was changed to that of Pd (from 195 to 106). 
There is a semi-quantitative agreement with the actual 57 Fe phonon partial DOS of 
Ll2-ordered Pd| Fe, suggesting the importance of the mass itself, but quantitative 
predictions of vibrational entropy are risky in this system, which undergoes an Invar 
transition under pressure. 

13. Unmixing 

13.1. Computational Studies of Unmixing 

The transformation from a disordered solid solution to two phases of distinct chemi- 
cal composition causes major changes in local chemistry, and substantial changes in 
configurational entropy. Nevertheless, there can be a large change in vibrational en- 
tropy when forming a disordered solid solution (Section 4.1). There are challenges in 
calculating phonons in disordered solid solutions. In first-principles electronic struc- 
ture studies by Ozolins, Wolverton and Zunger [209], effects of vibrational entropy 
were deduced (and effects on the miscibility gaps were assessed), by subtracting a 
calculated configurational entropy from an experimental total entropy of formation. 
This left a residual entropy of formation that was probably vibrational. For the 
formation of disordered equiatomic Ni-Au from pure Ni and pure Au, the vibrational 
entropy of formation was +0.48 fee /atom. For a solid solution of Cu-Ag, the vibra- 
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Figure 19. (a) Phonon partial DOS curves for 57 Fe, <7Fe(e), for Ll2-ordered Pt^Fe, measured by 
NRIXS and calculated using force constants of [205]. (b) Phonon partial DOS curves for 57 Fc, 
PFe(e), for Ll2-ordered Pd| 7 Fe, measured by NRIXS and calculated using force constants of [208]. 
(c) Phonon partial DOS curves for 57 Fe, gF e (e), for Ll2-ordered Pd| 7 Fe, calculated by using the 
force constants of Pt| 7 Fe [205] and the mass of Pd. After [207]. 



tional entropy of formation was +0.37 fee /atom. These are large contributions to the 
thermodynamics of unmixing. The authors suggested that the atomic size mismatch 
could be responsible for such large effects. From a computational study on unmixing 
in the Ni-Cr system, there was a suggestion that the discrepancy between phase 
boundaries between experiment and first-principles calculations without vibrational 
entropy might be useful for identifying the importance of vibrational entropy [210]. 

Zarkevich, Tan, and Johnson used first-principles methods to study the unmixing 
phase diagram in the alloy systems Ca-Sr, Au-Pt, Rh-Pd, Rh-Cu, Rh-Ag, and Rh- 
Au [211]. They also tested whether simpler methods could explain many of the 
trends, such as whether T cr i t could be predicted from a simple ratio of enthalpy to 
entropy. They discussed the difficulty of including vibrational entropy in simple 
estimates, but they showed some success with scaling the effects of vibrational 
entropy with differences in electronegativity of the elemental constituents. It may 
be difficult to use differences of electronegativity for quantitative estimates of the 
effects of vibrational entropy, especially for concentrated alloys, but it does seem 
that elements with large differences in electronegativity may show large effects of 
vibrational entropy in unmixing transformations. 
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Burton, van de Walle and Kattner performed first principles computational stud- 
ies of unmixing of three quasi-binary systems with the wurtzite structure, AIN-GaN, 
GaN-InN, and AIN-InN, using the bond-stiffness versus bond- length model [212]. 
With similar structures, there is a general expectation that the solubility should 
depend on differences in the ionic size of the metal atoms. The critical temperature 
for unmixing did correlate with the difference in ionic radii, but the effects of vibra- 
tional entropy did not. For example, although the phonon frequencies are expected 
to increase rapidly when larger In atoms are substituted on Al atom sites, the 
effects of vibrational entropy on the miscibility gap were surprisingly symmetrical 
in composition. This may be a consequence of the wurtzite structure, which is not 
close packed and may allow more local relaxations. 

13.2. Cluster Expansion Method for Solid Solutions 

A phonon DOS, a function of energy, pertains to a specific configuration of atoms. 
The interatomic forces are relatively short-ranged, so we expected that specifying 
local configurations of atoms would be sufficient to account for the phonon DOS. 
The cluster expansion method is a systematic way to construct a phonon DOS 
by summing contributions from different local clusters of atoms, starting with the 
chemical composition, and adding contributions from clusters of increasingly larger 
size. The cluster expansion method, which is extensible to DOS functions [213], is 
most useful when the phonon DOS can be accounted for with relatively few terms. 

A phonon DOS <7j(e) is related to the local atom arrangements in a struc- 
ture (the structure is designated j) through "interaction phonon DOS" functions, 
{Go(e), G\ (e), G2{e)...G v -\ (e)}, where v is the number of terms in the cluster ex- 
pansion. The ordered set of interaction DOS functions can be written as a vector, 
(5(e), so 

9j (e) = G(e) £ , (125) 

where £j is a set of atom correlation functions, {^}, in the specific composition or 
structure, j. With DOS curves from a set of different compositions or structures, 
g(e) = {gi(s), g2(s)...}, again written as a vector, we can solve for a set of different 
interaction DOS functions using an inversion method of Connolly and Williams 
[214], 

G(e)=g(e)C 1 , (126) 

where the matrix is obtained from the correlation functions below. 

The correlation functions £j2, • ■•} are obtained by geometrical analysis of 
atom positions in each structure, j. Specifically, for a cluster composed of n atoms: 

1 11 

^ = TrE^> ( 127 ) 
JV « i=i 

where pi are specific sites. Equation 127 is a product of "spin variables," where 
a = +1 for an A-atom and a = —I for a B-atom in a binary alloy. Consider, for 
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example, correlation functions for a disordered solid solution, where each site has 
the probability of xb of being occupied by a B-atom, and 1 — xb by an A-atom. The 
average for n = 1 is £i = <T£Xb+C/i(1 — xb) = 1 — 2xb- For clusters of two atoms (i.e., 
pairs), both sites are occupied randomly by A- or B-atoms, so £2 = (1 — 2xb) 2 , and 
for an n-atom cluster in a random solid solution £ n = (1 — 2xe) n . For an equiatomic 
alloy with xb = 0.5, £ n = for all n except n = 0. For n = 0, the interaction 
phonon DOS, Go(e), is the phonon DOS of the equiatomic solid solution, and all 
other G n (e) are corrections to this function for deviations from the equiatomic. 
The cluster functions G n (e) contribute to the DOS in proportion to how many n th 
clusters there are in a disordered alloy of composition c, compared to how many 
there are in an equiatomic disordered alloy. 

The same systematics can be used for the phonon partial DOS functions of A- 
and B-atoms. Here the index d denotes either "A" or "B." The analogs of Eqs. 125 
and 126 are 

9dj(e) = G d (e) £ dj (128) 
G d (e) = g d {e) if 1 (129) 

Neutron weighting (Sect. 8.3) is obtained with the phonon partial DOS (PDOS) 
functions 40 

ffr( £ )=E^^( M /^)^^), (130) 
d ma 

where each contribution to <7° w (e) (for the composition or structure j) is weighted by 
the phonon scattering efficiency of its A- and B-atoms and their partial DOS g^j{e)- 
The constant (M/a) provides normalization. It can be shown that the neutron- 
weighted interaction DOS functions can be written 

G™(e) = Y,RdC d U- 1 9d(e), (131) 
d 

where the matrix accounts for the neutron weighting [188, 215]. The matrix 
Rd relates the interaction DOS (IDOS) functions to the interaction partial DOS 
(IPDOS) function of species d (A or B) 

G{e) = Y J RdGd(e). (132) 
d 

The Rd is obtained from the correlation functions of the disordered solid solution. 

For a disordered solid solution, Lucas, et al. [188, 215, 216] found strategies to 
solve for the phonon partial DOS functions from neutron-weighted DOS curves of 



40 There are two subtleties: 1) It is assumed that the Debye- Waller factors for the two species of 
atoms do not differ significantly. This assumption should be adequate at low temperatures. 2) If 
the phonon spectra from all structures {j} do not span the same range of energy, the energy can 
be renormalized to the highest energy mode in each structure. This offers a good comparison of 
partial DOS functions, but requires further assessments for the entropy. 
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solid solutions. Arguably the best method is to obtain sets of experimental data 
with differences in the phonon scattering cross sections of the different species, 
such as with Fe alloys studied by neutron scattering and with 57 Fe inelastic nuclear 
resonant x-ray scattering (NRIXS). Alternative approaches can be developed with 
the assumption that the partial DOS of the host atoms is unaltered by the solute 
atoms, but this assumption can be violated if the solute atoms alter the electronic 
band structure or band filling. Improvements in the statistical quality of experimen- 
tal neutron-weighted phonon DOS functions offer other opportunities for comparing 
samples with rather similar chemical compositions, and further development of these 
methods is expected with improvements in neutron instrumentation. 

13.3. Fe-Cr Solid Solutions 

Using the cluster expansion method of the previous section, Lucas, et al., obtained 
a high degree of detail about how the vibrational entropy depends on local configu- 
rations of atoms in bcc Fe-Cr [188, 215]. The n = term, Go,o(£)) is the Cr partial 
DOS of a random solid solution of FesoCrsn. The n = 1 term gives the effect on 
the PDOS of having different numbers of Fe and Cr atoms in the nearest-neighbor 
(NN) shell of the Cr atom. 41 The function Gici{s) is added to Go,cr(£) when the 
atom in the NN shell of the Cr atom is Cr, and subtracted when it is Fe. Figure 20 
shows that the effect of Cr atoms in the NN shell of a Cr atom is to decrease the 
number of modes in the ranges of < e < 0.59 and 0.83 < e < 1, while increasing 
the number of modes in the range 0.59 < e < 0.83. The triangle term n=2 gives 
the effect of the like or unlike pairs of atoms within the NN shell of a Cr atom. The 
function G2,Cr(e) adds to both Go,Cr(e) and G^ot^) when a pair in the NN shell 
of the Cr atom is composed of two of the same atoms. Pairs of like atoms decrease 
the number of modes in the ranges < e < 0.56, 0.80 < e < 0.88, and 0.94 < e < 1 
while increasing the number of modes in the remaining ranges, whereas unlike pairs 
do the opposite. 

The IPDOS functions can be used to construct the PDOS curves of Fe and Cr 
atoms in alloys (Fig. 21). These PDOS curves can then be used to construct the 
neutron-weight-corrected DOS curves. The neutron-weight correction stiffens all of 
the DOS curves because the Fe PDOS curves have more low-energy modes than 
the Cr PDOS curves, especially below e = 0.7, but Fe atoms are more efficient 
scatterers. For Fe-Cr alloys, the PDOS curves of both Fe and Cr atoms undergo an 
average softening upon alloying, giving a net positive vibrational entropy of mixing. 
The softening of the Cr PDOS curves is larger, so at low concentrations of Fe in Cr, 
the phonon entropy of mixing rises relatively rapidly with Fe concentration, giving 
a skewness to the concentration dependence of the vibrational entropy of mixing. 

Cluster expansion methods were originally developed for thermodynamic quan- 
tities such as entropy, and the same notation is useful for interaction entropies of 
solid solutions, s"p e and s" Cr - For a random solid solution, the vibrational entropy 
of mixing is, in terms of the interaction phonon partial entropies 

A<T ib = -4X0,(1 - X Cr ){[sf h Cr + <Fc]/2 + *Cr sf h Cl + (1 - ZCr^Fe} ■ (133) 

41 Lucas, et al., used a combined first and second nearest-neighbor shell of the bcc structure as 
the NN shell [188, 215]. 
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Figure 20. Interaction phonon DOS functions of (a) Fe and (b) Cr in bec solid solutions of Fe-Cr 
alloys, determined from neutron scattering (solid lines), and direct measurement by 57 Fe NRIXS 
(points). The n = 1 and n = 2 curves are offset vertically. Note that the energy scale was 
renormalized to the energy of the highest energy mode in each structure. 



The function —Axq t (1 — xq t ) is a parabola, scaled by the average [s^ r + ^i'fcV^- 
This average is negative, so the term introduces a positive entropy of mixing, with 
a maximum at xo = 0.5. The remaining two terms associated with are the 
equation of a line, giving an asymmetry to the vibrational entropy of mixing when 
s 2Ci 7^ ^sTfc The vibrational entropy of mixing is shifted toward pure Cr because 
G ! 2,Cr(£) < G2,Fc(£) in the energy range to 25 meV. The maximum critical tem- 
perature and composition of the miscibility gap in Fe-Cr are 905 K and xc r =0.51, 
respectively. Assuming that the values of sj^ for n > do not change appreciably 
with temperature, Eq. 133 can be used to determine the effect of the vibrational 
entropy on the miscibility gap for spinodal decomposition. When the vibrational 
entropy of mixing in the bec phase is removed, the miscibility gap for a random 
solution of Fe-Cr is shifted up in temperature and toward pure Cr, with a critical 
temperature and composition of 1268 K and xo = 0.61, respectively [188, 215]. It is 
known that the enthalpy of mixing is also shifted toward pure Cr [217, 218, 219, 220], 
so the vibrational entropy counteracts the asymmetry of the enthalpy of mixing, 
especially at higher temperatures. In short, the vibrational entropy of mixing lowers 
the miscibility gap and shifts its maximum towards the equiatomic composition. 
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Figure 21. PDOS curves from the cluster expansion fitting method for the concentrations of the 
alloys measured: (a) Fe PDOS curves (solid lines), and measured pure Fe DOS (dashed line), and 
(b) Cr PDOS curves (solid lines), and measured pure Cr DOS (dashed line). After [188, 215]. 

14. Interstitial Alloys 

14.1. Interstitial Oxygen and Nitrogen 

Inserting small atoms into interstitial sites of a crystalline host matrix is sometimes 
treated as a simple problem of distributing the small atoms over a set of independent 
lattice sites. The energy of insertion of each atom is assumed a constant, and a config- 
urational entropy can then be calculated readily with a form of Eq. 11, for example. 
The real situation is more complicated, owing to distortions of the matrix around 
each interstitial atom. Sometimes the elastic interactions between distorted unit cells 
cause structural organizations of the interstitial atoms. Anisotropic structures such 
as martensite in steel are formed this way, accompanied by large, long-range elastic 
distortions. These distortions are expected to alter the vibrational frequencies of 
both host and interstitial atoms, and the vibrational entropy. Such effects are not 
yet well understood. 

Danilkin, et al., performed inelastic neutron scattering studies on V, Nb, and Ta 
with concentrations of up to llat.% of interstitial O and N [221, 222, 223]. They 
argued that the mass difference between the host and interstitial elements causes 
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separation of the phonon DOS into bands dominated by the partial DOS of the metal 
atoms at low energies, and dominated by the partial DOS of the interstitial atoms 
at high energies. Although a full justification is not presented, it may be reasonable 
to neglect the contribution from the interstitial partial DOS in the low-energy part 
of the spectrum, and ascribe all changes to the phonon partial DOS of the metal 
atoms. With this assumption, large changes in vibrational entropy of the host metal 
are obtained. For example, in a V-3at.%0 alloy, the average vibrational frequency 
of the V atoms increases by approximately 6%, which is a change of +0.18 &b/V 
atom with respect to pure V. This is large compared to the configurational entropy 
of mixing for an alloy of 3at.% O, which from Eq. 11 is — 0.135 k^/O atom. It is also 
interesting that these effects of interstitial atoms on the vibrations of host atoms 
are sub-linear at higher solute concentrations. Similar, although somewhat smaller, 
effects were found for nitrogen in vanadium, and for oxygen in niobium. For alloys 
of Ta-N, for which the mass difference is quite large, the effects were somewhat 
smaller, although the authors point out their concern that the measurements did 
not average the coherent scattering from Ta over a sufficiently wide range of Q. 

14.2. Hydrogen and Hydrides 

Hydrogen and deuterium differ in mass by a factor of two. The characteristic vibra- 
tional frequencies of H2 and D2 molecules differ by almost exactly y/2. Analogous 
mass effects are expected for hydrogen in a metal hydride. Owing to its low mass, 
hydrogen typically has high vibrational frequencies, which can show effects of quan- 
tized occupancy at temperatures well above room temperature because the Planck 
occupancy factor is still small. The V-H and the V-D phase diagrams have been 
measured, and the phase boundary between the disordered solid solution phase and 
the two-phase region is shifted to lower temperature for V-D by approximately 15% 
[224]. 42 With deuterium there is a greater difference in the vibrational entropy of 
the solid solution compared to the unmixed state at lower temperatures, and this 
difference with respect to hydrogen seems understandable as a mass effect [225]. 
It has also been noted that with substitution of deuterium or tritium, the lattice 
parameter of LiH is reduced at low temperatures, although this difference decreases 
with temperature [226]. At high temperatures in the classical limit, it is expected 
that the mass effects on the thermal expansion of LiH versus LIT will disappear. 

The alkali hydrides LiH and NaH both have the Bl (NaCl) structure at ambient 
pressure, but NaH transforms to the B2 (CsCl) structure at a pressure of 29 GPa, 
and LiH is predicted to transform to B2 at 308 GPa [227]. A computational study 
of phonons in Bl LiH and LiD by Yu, et al. [227], showed a large isotope effect 
on the three high-frequency optical branches that have a large partial DOS for 
H atom motions. A much smaller isotope effect was found for the three acoustic 
branches. The calculated heat capacity of Bl LiH was in good agreement with exper- 
imental measurements, which also reported that the Griineisen parameters increase 
at lower temperatures [228]. Since the phonon spectrum extends to approximately 
100 meV for LiH, the difference in the zero point energies of LiD and LiH at ambient 
temperature is important for the vibrational thermodynamics, and there are large 



There are also differences in the phase fields of the V-H and V-D phase diagrams. 
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zero-point effects on the lattice parameters and bulk moduli. The calculated phonon 
frequencies of LiH and NaH varied with volume, and these frequencies were used with 
the harmonic entropy expressions (i.e., a quasiharmonic approximation was used) 
to calculate the Helmholtz free energy. The calculated Bl to B2 transformation 
pressure for NaH was in good agreement with experimental measurement [227]. 

Equilibrium between a system containing a hydrogen sorbent in contact with a 
gas reservoir at fixed temperature and pressure ensures that the enthalpy of sorption, 
AH, is dissipated as heat. 

AH = H g -H s = T(S g - S s ) . (134) 

with the subscript g for the full system with hydrogen in the gas, and s for the same 
system with hydrogen in or on the solid sorbent. The entropy of hydrogen gas is 
large, 5 gas = 131 J(mol K) -1 = 7.9/eB/atom at STP. For the same pressure of gas, 
sorbents with larger AH tend to have higher operating temperatures. 

A consequence of Eq. 134 is the release of heat during hydrogen sorption. Consider 
the case of lithium borohydride, for which AH = 67kJ/mole at 400° C. The heat 
release from the sorption of 2500 moles of hydrogen is 167 MJ. For a refueling in 3 
minutes, one megawatt of heat must be dissipated! Reducing this heat load requires 
either lower sorption temperatures, or a reduction of S g — S s . The degrees of freedom 
of diatomic hydrogen gas that contribute to S g are translational, rotational, vibra- 
tional, electronic, and nuclear. The first two contributions are big, and generally 
unavailable to hydrogen in a solid, reducing the maximum value of S s . On the other 
hand, the vibrational frequency of the H2 molecule has a characteristic temperature 
of 6215 K, suppressing its contribution to the entropy at ambient temperature. It is 
possible to get lower vibrational frequencies in compounds. Ozolins, Majzoub, and 
Wolverton [229] found effects of this sort, and suggested an approach for minimizing 
S g — S s . Storing hydrogen in complex anions such as NH2 and BH4 can produce 
crystal structures with low- frequency rotational modes, which increase the vibra- 
tional entropy of S s . After the hydrogen is removed, it is helpful if the solid has a 
low vibrational entropy owing to stiff phonons, further reducing S g . Nevertheless, 
with S gas = 7.9 fee/atom, it is a challenge to decrease S g — S s substantially. 

14.3. Lithium Intercalation and Insertion 

There are numerous electrochemical studies of the insertion and de-insertion of 
lithium in materials. A few of these studies have reported precise information on 
the chemical potential of the Li atoms, and how it changes with concentration. In 
equilibrium (i.e., steady-state without changing the charge of the electrochemical 
cell) , the difference in voltage of the anode and cathode of a reliable electrochemical 
cell is precisely the difference in chemical potential of the Li atom between one 
electrode and the other. A counter-electrode of pure Li metal is especially useful 
for thermodynamic measurements, because there is no difference in the chemical 
potential for the different Li atoms as they are deposited on, or removed from, 
a pure Li metal electrode. By varying temperature, it is possible to measure the 
temperature dependence of the free energy per atom, and hence determine the 



Vibrational Thermodynamics of Materials 



71 



entropy of intercalating additional Li into the material 

AS(x) = F^ , (135) 

where F is the Faraday constant and Eq is the equilibrium voltage of the cell. 
Although Eq may be volts, and may change with temperature by only 0.1 mV per de- 
gree K, accurate measurements of AS(x) are practical with modern instrumentation 
and good electrochemical cells. 

14.3.1. Lithium Intercalation in Graphite 

Graphitic carbons are in widespread use as anode materials in rechargeable Li 
batteries. Under ambient conditions, Li atoms can be inserted into graphite to a 
concentration of LiiCg. Several intermediate phases are formed. As the concentration 
of Li, x, increases in Li x Cq, the first phase is a disordered distribution of Li atoms in 
the van der Waals gap between graphene planes. Figure 22 shows a large and positive 
entropy of insertion for small x. Assuming that all sites in the disordered phase are 
filled at x = 0.25, the concentration dependence of Eq. 11 predicts the curved part of 
AS^ic which is in reasonable agreement with the experimental results for A5 meas . 
Over this range of composition, AH meas is reasonably constant, so perhaps a lattice 
gas model can successfully explain the intercalation thermodynamics. 

It is known that further lithiation occurs in stages. A Stage II compound (with 
Li atoms filling alternate gaps between graphene planes) is formed up to the con- 
centration x = 0.5, followed by a Stage I compound (with Li atoms filling all 
gaps between graphene planes). These compounds form by first-order reactions, 
with an approximately constant enthalpy and entropy for each additional amount 
of transformed material. Figure 22 does show that AS mcas and AH meas are fairly 
constant, with a break at the boundary of x = 0.5 that marks the onset of Stage I 
compound formation. 

Figure 22 shows that the entropies of formation of Stage II and Stage I compounds 
are negative. Experimental studies of the vibrational spectra of intercalated lithium 
atoms provide two Debye temperatures, a #d_l — 893 K for the one dimension of 
lithium vibration normal to the graphene planes, and a D y ~ 392 K for the two 
dimensions of lithium vibrations in the plane of the graphene sheets [230] . With a 
Debye temperature 6>do — 380 K for bec Li metal [231], the vibrational entropy of 
lithiation A5 ca i c is obtained as the flat line in Fig. 22. In essence, the vibrational en- 
tropy of intercalation is negative because the modes of Li vibration perpendicular to 
the graphene planes are highly constrained, and have high frequencies. The resulting 
negative vibrational entropy of intercalation dominates the entropy. There should 
be little change in configurational entropy when Stage I or Stage II compounds are 
formed because they are ordered. 

14.3.2. Lithium Insertion in Li x CoC>2 

The layered hexagonal compound LLjCoO^ preserves a transition-metal oxide crystal 
framework as Li + ions are inserted or extracted. Because the chemical potential of 
the Li atom is lower by approximately 4 eV in Li x Co02 than in graphite, these oxides 
are excellent for service in rechargeable Li batteries. The "03" hexagonal structure 
is stable over the composition range, ~ 0.58 < x < 1, the approximate range used 
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Figure 22. Composition-dependence of the entropy and enthalpy of intercalation, AS(i) and 
AH(x). Labels ASmcas and Aff mcas refer to the experimental data. The label AS'caic refers to 
the calculation of the entropy of intercalation (see text). After [232]. 



in practice for rechargeable batteries. For 0.95 < x < 1, however, the 03 phase has 
a lattice parameter that is smaller than for x = 0.833, and there is a region of two- 
phase equilibrium between x = 0.833 and x = 0.95. There is a concomitant change 
from metal to insulator between x = 0.833 and x = 0.95. In the composition range 
0.58 < x < 0.833 there is a disordered distribution of Li atoms on the available sites, 
so the configurational entropy is expected to vary strongly with composition. 

Phonon DOS measurements were performed on samples of Li x Co02 for various 
values of x in the 03 structure, including x = 0.62 and x = 1.0 [233]. Perhaps 
surprisingly, there was little difference in the DOS or the vibrational entropy of 
these two compositions. On the other hand, there is a considerable difference in the 
vibrational entropy of the Li in Li metal versus in Li x Co02- Because the phonon 
frequencies of Li metal are much lower than those of Li x Co02, as Li is transferred 
from metal to oxide, the vibrational entropy decreases. The contributions to the 
entropy of lithiation of Lia;Co02 are shown in Fig. 23. The data points were obtained 
with measurements that exploited Eq. 135. The smooth curve labeled A5 c f was 
calculated from the configurational entropy of the disordered solid solution (it is the 
derivative dS c f/dx because the open circuit voltage measures the change in entropy 
as additional lithium is inserted in the material). In Fig. 23, A5 c f was lowered 
vertically to make a best fit to the electrochemical results. This offset is mostly 
accounted for by the vibrational entropy of lithiation, ASVib, which is large and 
negative, approximately —2.4 fee /atom, although it has only a small variation with 
x. (The correction from the electronic entropy, AS e \, is not large because there are 
few states near the Fermi level in these materials.) In summary, and in analogy to Li 
intercalation into graphite, the Li + inserted into Li x Co02 is more stiffly constrained 
than in Li metal, giving higher vibrational frequencies and a lower vibrational en- 
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Figure 23. Configurational entropy AS c f (offset vertically), vibrational entropy AS'ph, (and elec- 
tronic entropy AS c i correction) of the lithiation reaction in LLcCoCh over the composition range 
0.6 < x < 0.833, superimposed on electrochemical data. 



tropy. Since this does not change much with Li concentration, the dominant effects 
probably involve the Li and the host structure (rather than Li-Li interactions). 



15. Martensite and Twinning 

Iron-carbon martensite is a classic interstitial solid solution. The deformation po- 
tential for C in bet Fe is of order 1 eV - much higher than for most substitutional 
solid solutions [234]. With a carbon concentration of 4at.%, the tetragonality of the 
martensite unit cell is nearly 4%, which causes a large distortion of the nearest- 
neighbor bonds between Fe atoms and between Fe and C [235, 236, 237]. Effects 
on vibrational entropy from rearrangements of carbon atoms are plausible. The 
phonon DOS might be changed substantially by the presence of C in Fe, perhaps 
with a stiffening of the host Fe DOS as was found in the work by Danilkin, et 
al., for V, Nb, and Ta with interstitial atoms (Section 14.1). On the other hand, 
since C diffuses freely in bec Fe, we should expect a substantially weaker bonding 
between Fe and C than between transition metals and O, reducing any effects of 
stiffening. The thermodynamics of martensitic transformations has a long history 
of study [238, 239, 240, 241, 242, 243, 244]. Nevertheless, most of the interest has 
been on the thermodynamics of the martensitic transformation itself, as discussed 
below, and not of the thermodynamics of formation of the interstitial solid solutions 
of Fe-C. 

A structurally stable crystal must provide restoring forces against shear and 
dilatation strains that would change its shape. The stability of a cubic crystal 
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requires these conditions on the elastic constants 

C n > , C 44 > , Cfi > C\ 2 , i(Cn - C12) > . (136) 

Each of these stability criteria corresponds to a different mode of deformation. For 
example, the last inequality, of particular interest in the study of martensite, defines 
the stability with respect to small amplitude shears on (llO)-type planes in (110) 
directions. If the ratio (Cn — C12) vanishes at a given temperature, the crystal 
will shear spontaneously in a "soft mode instability." The earliest explanations of 
martensitic transformations [114, 245] were based on soft mode instabilities. Clapp's 
[246] formulation used the the free energy of a crystalline solid 

-^6 6 -^666 

F = F + - djViVj + ^jEEE C ijk THVjTlk + - (137) 

i=l j=l ' i=l j=l k=l 

where Fo is the free energy of the strain-free material, and Voigt notation for strains 
rj is used. 43 Stability with respect to strain is governed by the sign of each of the six 
eigenvalues of the 6x6 matrix of second derivatives: 

iB^r F " w_1,8) - (138) 

If a particular eigenvalue vanishes at a given temperature, the material is unstable 
with respect to the deformation mode of the associated eigenvector. The stability 
analysis applied to a free energy written in terms of only the first summation in Eq. 
137 will reproduce the inequalities presented in Eq. 136. 

It is a general feature of the bcc structure that it is weak against longitudinal 
waves with 2/3(111} wavevectors. Phonon dispersion curves of all bcc elements 
show a dip in frequency of the transverse dispersion near this wavevector, or at 
least an anomaly. (In contrast, phonon dispersions of fee crystals tend to increase 
monotonically to the edge of the Brillouin zone.) Many bcc crystals have phonon 
dispersions that also soften at specific frequencies along the transverse [110] branch, 
but the specific wavevector of these soft modes varies for different alloys. At temper- 
atures near some martensitic transformations, these low-frequency transverse modes 
dip towards zero frequency, becoming unstable. This behavior could originate with 
a temperature-dependent change of second nearest-neighbor force constants [47], 
because the bcc structure is not stable with first-neighbor central forces alone. With 
large atom displacements, however, the interatomic force constants are not expected 
to behave as for small displacements, and the mechanics of how a soft mode can 
lead to structural changes is not yet well understood. 

Many bcc metal crystals are unstable energetically at low temperatures, where 
a close-packed structure is preferred. In pure elements, the bcc-fcc transformation 
cannot cause any change in configurational entropy because the shear and shuffle 
processes of atom movements create the product phase deterministically. Owing 
to the weakness of bcc structures against various distortions, Friedel suggested 
some time ago that soft transverse modes stabilize the bcc structure at elevated 



43 Here rji =en,= £22,^3 = £33, V* = £ i2,?75 = £13 and ije = e 2 3- 
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temperatures [247] (also discussed in Section 10.1). Other considerations of entropy- 
changes during martensitic transformations and twinning in alloys were described 
by Cahn [248]. 

Murakami studied the cubic to orthorhombic martensitic transformation in Au- 
Cu-Zn alloys [249]. Elastic constants were obtained by ultrasonic measurements, and 
Debye temperatures were estimated from them. The importance of vibrational en- 
tropy was demonstrated, at least qualitatively. Tong and Wayman presented similar 
ideas about the fee to hep martensitic transformation in Co-Ni [242]. 

Planes, Ortm, Mahosa, et al., performed a number of calorimetric studies of the 
heats of martensitic transformations in shape-memory alloys, which tend to have 
good thermodynamic reversibility. From the measured heats from Cu-Al-Zn, Cu-Al- 
Ni, and Cu-Al-Be alloys, AH, and the temperatures of transformation, the entropy 
was obtained as AS = AH/T. A typical value of AS was 0.2 k-Q /atom [250, 251, 
252]. The size of AS increased with the electron/atom ratio of the alloys, with AS 
varying by about 15% over the composition ranges studied. These investigators also 
showed that when the transformation temperature is low with respect to the Debye 
temperature of the alloy, the vibrational entropy of the transformation is suppressed 
- a consequence expected from a low occupancy of most phonon modes. Finally, they 
presented evidence that the martensitic transformation temperature depends on the 
presence of low-energy modes in the high temperature phase that help stabilize it, 
and the transformation is promoted when these modes are depopulated at lower 
temperatures. 

The martensite-to-austenite transformation in NiTi has a latent heat of (0.5 ± 
0.05) fee/atom, as measured by calorimetry [253]. The phonon DOS of the two phases 
were obtained from inelastic scattering spectra after using a lattice dynamics simula- 
tion to correct for the different phonon scattering weights of Ni and Ti. These phonon 
DOS curves of the austenite and martensite, Fig. 24, show substantial differences. 
From these data, it seems that all the transformation entropy from calorimetry 
is vibrational entropy. From simulations and the measured DOS, the vibrational 
entropy of austenite was found to be larger than that of martensite because the 
transverse acoustic (TA(2)) and longitudinal acoustic modes of austenite are softer 
than those of martensite. (Similar changes in acoustic modes of NiTi were reported 
previously [254].) Simulations suggest that this originates with a large and negative 
INN transverse force constant, which may also be responsible for the soft modes 
involved in the mechanism of the martensitic transformation. 

The situation for martensitic transformations in Fe alloys is more complicated. 
Substantial changes in magnetic entropy typically accompany the martensitic trans- 
formation. For Fe-Ni alloys, the martensitic transformation is from fee at high 
temperatures to approximately bec at low temperatures. Figure 25 shows phonon 
DOS curves from Fe-29 at.% Ni [255]. The phonon DOS curves were corrected for 
neutron weighting, using the partial DOS of 57 Fe atoms measured by Franz, et 
al., [256]. The vibrational entropy change of the martensitic transformation was 
ASVib = 0.12&B/atom, with the fee phase having the larger vibrational entropy. It 
is interesting to see that the Ni atoms have relatively low vibrational frequencies in 
the bec phase, reducing the vibrational entropy of the transformation. This ASVib 
was inconsistent with the latent heats of the forward and reverse martensitic trans- 
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formations, which were not equal. The origin of the discrepancy was the change in 
magnetic entropy. The martensitic transformation is hysteretic, and occurs below 
the Curie temperature of the bcc phase for the forward fcc^bcc transformation, 
and above this temperature for the reverse transformation. The change in magnetic 
entropy of the martensitic transformation in Fe 29 at.% Ni proved to be about twice 
as large as the change in vibrational entropy. 

The Heusler alloy N^MnGa exhibits a "magnetic shape memory" effect, where 
a magnetic field can change the shape of the material much as temperature changes 
the shape of conventional shape memory alloys like Ni-Ti. The alloy Ni2MnGa also 
changes its shape by a reversible martensitic transformation, and there have been 
several studies of phonons in the cubic "austenite" and tetragonal "martensite" 
phases. Several soft phonon anomalies near the transition temperature were observed 
in the phonon dispersion curves, and it was suggested that they originate in part 
from electron-phonon interactions of a nested Fermi surface [257]. Further work 
on Ni2MnGa showed that the transition temperature depends on the electronic 
structure [258], and the thermal excitations of both magnons and phonons [259]. 
The magnetic and vibrational effects on the martensitic transition in Ni2MnGa 
seem to be comparable in importance. 



16. Microstructural Contributions 



The entropy is typically largest when the number f2 of Eq. 1 is obtained by counting 
entities with about one degree of freedom per atom. Consolidating n degrees of 
freedom into one reduces S by a factor of n, so mesoscopic and microscopic structures 
in materials generally have too few degrees of freedom to compete with entropic 
contributions from atomistic quantities. Usually it is the energy that depends on 
the microstructural degrees of freedom, but some exceptions are presented here. 
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Figure 25. (a) Partial phonon DOS curves for Fe and Ni in the martensite, at room tempera- 
ture. Line: Fe-pDOS obtained by Franz, et al. [256]. Dots: Ni-pDOS determined from analysis of 
data of [255, 256]. (b) Room-temperature phonon DOS curves. Dashed line: neutron- weighted 
DOS for martensite. Diamonds: neutron-weight-corrected DOS for martensite. White circles: 
neutron- weighted DOS for austenite (the neutron- weight in this case is minor). 

16.1. Elastic Energy 

If the thermal expansion of a crystal is isotropic, all distances will expand by the 
same fractional length. Even an irregular packing of crystals in a polycrystalline 
microstructure will expand uniformly, and will not develop internal stresses with 
changes in temperature. On the other hand, when thermal expansion is crystallo- 
graphically anisotropic, an irregular polycrystalline microstructure develops internal 
stresses as adjacent crystals, differing in crystallographic orientation, try to ex- 
pand by different amounts along different directions. 44 Here we consider the case 
where a stress-free polycrystalline microstructure is prepared by annealing at high 
temperature, but anisotropics of thermal contraction create internal stresses upon 
cooling. 

The DO22 structure is similar to a cubic fee structure, but has a compositional 
modulation along (001) that gives a tetragonal distortion to the unit cell. If the 
DO22 structure is formed in equilibrium at high temperature and cooled, a polycrys- 
talline microstructure will develop internal stresses at low temperature. Assuming 
the strains are all elastic, these internal stresses are relieved as the microstructure is 
heated to modest temperatures, and this stored elastic energy is converted to heat. 



Some ceramic materials with unit cells of low symmetry may fracture under these internal 
stresses upon cooling. 
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Figure 26. Dots: differential heat capacity of two samples of N13V. Sign is partially-disordered fee 
minus D0 22 -ordered. After [260, 196]. 

It therefore takes less heat input to raise the temperature of this microstructure, 
and its heat capacity is suppressed. Evidence of this suppressed heat capacity was 
found for DO22 N13V and PdsV alloys [260, 261]. In these experiments, differential 
scanning calorimetry was used to compare alloys in two phases: 1) cubic, prepared 
by quenching the sample from a high temperature in an fee region of the phase 
diagram, and 2) DO22 ordered, prepared by annealing in the DO22 region of the 
phase diagram. Figure 26 shows the differential heat capacity of these two samples 
- the partially disordered cubic sample had the larger heat capacity. The shape 
of the differential heat capacity curve is quite different from the shape of Fig. 8b, 
however. The curve labeled "Microstructure" in Fig. 26 is the contribution obtained 
by assuming that the DO22 ordered sample retains internal stresses at low temper- 
ature that are relieved as the solid expands. The curve labeled "DOS" in Fig. 26 
was obtained from experimental phonon spectra measured for the two samples by 
inelastic neutron scattering. The sum of the two contributions better accounts for 
the measured differential heat capacity than either contribution alone. 

Nevertheless, the fitted shape in Fig. 26 leaves much to be desired. The re- 
sponse to thermal stresses is probably more complicated than linear elastic behavior. 
Over a temperature range of 100 K, differences in thermal expansion along different 
crystallographic directions of 10 _5 /K generate strains in excess of 0.001. Over a 
temperature range of a few hundred K, local plastic deformation can be induced by 
anisotropic thermal strains, especially at locations of stress concentrations. Plastic 
deformation requires mechanical work to be performed by thermal expansion, in- 
creasing the apparent heat capacity (as opposed to decreasing it by relief of internal 
stress). The microstructural contribution to the heat capacity can therefore change 
sign, and plastic deformation can cause hysteretic behavior. 

The effects of elastic and plastic deformation during thermal contraction and 
expansion were studied in polycrystalline uranium, and its differential heat capacity 
was measured against a single crystal of uranium [262]. Both samples are expected 
to have nearly the same phonon entropy, so the study was able to isolate the mi- 
crostructural effects on the heat capacity. Powder diffractometry was used to obtain 
internal strain distributions along different crystallographic directions. These diffrac- 
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tion results were used to calculate the elastic energy in the microstructure, and this 
analysis of diffraction results was in reasonable agreement with the measured heat 
capacity. The differential heat capacity between polycrystalline and single crystal 
uranium changed sign with temperature as the relief of elastic stress upon heating 
gave way to plastic deformation. 

These effects are probably not best described as part of the intrinsic thermody- 
namics of crystalline phases, but rather a part of the mechanics of the microstruc- 
ture. They do affect bulk thermodynamic properties measurements, however, and 
such measurements offer a means for studying them. 

16.2. Nanostructures 

f6.2.1. Energy and Entropy of Nanostructures 

Nanostructured materials are of widespread interest in many fields of science and 
engineering. A practical question is whether nanostructures are sufficiently stable 
at their temperatures of service, and a basic question is how the thermodynamics 
of nanostructured materials differs from conventional bulk materials. In short, the 
surfaces or interfaces of nanostructures raise the energy. On the other hand, the free 
energy of nanostructured materials could be reduced by sources of entropy associated 
with configurations of atoms or excitations in the structure, helping to stabilize the 
material at finite temperatures. Consider, for example, a nanostructured material 
with a grain boundary enthalpy of 3kJ/mole [263, 264]. At a temperature of 500 K, 
an entropy of 0.75kB/atom (6.2 J(molK) -1 ) would null the free energy difference 
between the nanostructured material and the bulk material. 

The extra configurational entropy associated with a nanocrystalline microstruc- 
ture is probably not large. For consolidated nanocrystals, there can be a substantial 
fraction of atoms in grain boundaries, which are typically of 1 nm in width [264, 265]. 
Atoms in the grain boundary accommodate the mismatch of the bounding crystals, 
and cannot be positioned as arbitrarily as in a glass. Using a glass transition entropy 
of ~ Ike/atom would be a generous upper bound. This is shared over all atoms in 
the solid, however, so it should be reduced by the grain boundary fraction of atoms. 
The positions and orientations of the nanocrystals give other configurational degrees 
of freedom. Considering each nanocrystal as having nucleated at a point and grown 
to impingement with other nanocrystals, the configurations can be counted as the 
number of ways of distributing the nucleation points over a lattice of approximately 
atomic spacing. The orientational degrees of freedom can also be counted as rotations 
of f-atom increments at the edge of the nanocrystal. These configurational degrees 
of freedom are not large for a nanocrystal of 7nm dimension and 10 5 atoms, when 
averaged over all atoms in the material. (There will be additional degrees of freedom 
if the positions of grain boundaries are decoupled from the centers of the nanocrys- 
tals.) Perhaps 0.1 ke/atom can be obtained for the excess configurational entropy of 
such nanocrystals, but this seems generous. In summary, the configurational degrees 
of freedom of atoms become consolidated into collections of 10 atoms, reducing the 
entropy by this large factor, but there are nanostructural features that bring some 
new degrees of freedom. 

The vibrational entropy of nanostructured materials is less well understood. 
From the measurements available to date, the vibrational entropy difference be- 
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tween nanostructured and bulk materials is larger than the configurational entropy- 
difference, but it seems half as large as necessary to stabilize thermodynamically a 
nanocrystalline material at a temperature of 500 K. Nevertheless, vibrational entropy 
makes a difference for the stability of nanocrystalline materials. Furthermore, the 
vibrational dynamics of nanostructures are of interest in their own right. 

16.2.2. Low-Energy Vibrations of Nanostructures 

The eigenmodes of vibrations of small crystals have been calculated both for isolated 
crystals [266, 267, 268]. and for crystals with rigid constraints at their surfaces 
[269]. A distinguishing feature of these calculations of the DOS, later borne out 
by experiments, is an enhancement in the number of modes at low energies. A 
continuum analog for some of these low-energy modes are Rayleigh modes over 
the surface, which become more distinct from bulk atom motions with increasing 
mismatch of the acoustic impedances of the crystal and the medium in which it is 
embedded. These modes involve surface atoms, of course, but they also contribute to 
the partial DOS of atoms beneath the surface. A computer simulation of a periodic 
array of small nanocrystals with Lennard-Jones interatomic potentials also showed 
an enhancement of the low-energy modes of the phonon DOS [270]. The calculated 
DOS was similar to that of a glass, but curiously, the glass was more stable than 
the nanocrystalline microstructure of this simulation. 

Experimental measurements have shown a number of differences in the phonon 
DOS of nanocrystalline and bulk materials. In general, the energies of the van Hove 
singularities are the same for materials with nanocrystals or with large crystals, 
indicating that the interatomic forces are the same for most atoms situated locally 
in regions of good crystal [271, 272, 273, 274]. This is seen in the phonon DOS curves 
of Fig. 27 - although the DOS is broadened with smaller crystallite size, the energies 
of the van Hove singularities do not shift substantially [275]. This remains true as 
the nanocrystalline material is subjected to pressure, as shown for bcc Fe in Fig. 
28 [276]. At pressures of 0.2 and 9.5 GPa both bulk and nanocrystalline samples 
were bcc, but both had transformed to hep Fe at the highest pressure of 25-27 GPa. 
The energies of the van Hove singularities changed with crystal structure and with 
pressure, but these changes were similar for the nanocrystalline and bulk materials. 

It is clear from both Figs. 27 and 28 that the phonon DOS in the nanocrystalline 
material has more spectral weight at energies below about 20meV, and further 
analysis of the data of Fig. 27 showed that this enhancement at low energies scales 
inversely with the crystallite size. This enhancement of the DOS at low energies 
has been found in all studies of compacted nanocrystals [272, 273, 274, 275, 277, 
278]. This trend implies that the extra spectral weight at low energies originates 
from degrees of freedom in the microstructure, which may include vibrations that 
tend to localize at grain boundary regions or propagate along them, motions of the 
nanocrystals with respect to each other, or some combinations of these degrees of 
freedom. 

For consolidated nanocrystals, the enhancement in the phonon DOS at low en- 
ergies does not extend to arbitrarily low energies and aritrarily long wavelengths. 
In measurements of phonons in the fieV energy range [279], it was found that the 
enhancement in the phonon DOS around 10 /ieV was only about half of what was 
found around lOmeV. Some enhancement of the phonon DOS is expected owing 
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to lower average elastic constants in any composite material with elastically "soft" 
atom configurations at grain boundaries, and this may explain the suppression of 
the DOS in the /jeV energy range and the reduction in bulk elastic constants. The 
behavior at ^eV energies is therefore understood as a continuum limit, so the extra 
enhancement in the meV energy range is more interesting. It is likely caused by 
the extra degrees of freedom of the nanocrystalline microstructure. These probably 
originate with the spatial dimensions of the nanocrystalline features, of order 10 nm, 
and these effects therefore do not extend to the continuum limit. We note in passing 
that this behavior of an enhancement in the DOS at low, but not arbitrarily low, 
energies is reminiscent of the "boson peak," of glassy materials, hinting that the 
boson peak may originate from nanostructured features of the glassy structure. 

The enhancement of spectral weight at low energies is a widespread phenomenon, 
and is larger when the nanocrystalline dimensions are smaller. Further insight can 
be obtained by assessing whether these low energy excitations show characteristics 
of reduced dimensionality. In the continuum limit, all sound waves have linear dis- 
persions with e/h = lo oc k, so in three dimensions g(e) oc e 2 . Especially for DOS 
curves measured by inelastic nuclear resonant scattering, which offers the advantages 
of being both incoherent and with a low background, the exponent n of g{e) oc e n 
can be obtained by fitting g(e) as in Fig. 28. Doing so for three DOS curves in 
the bec phase at different pressures gave n = 2.00(±0.015) for the bulk material, 
and n = 1.95(±0.015) for the nanocrystalline material. Some excitations in the 
nanocrystalline material may have 2-dimensional character. The nanocrystals used 
for the data of Fig. 28 had grain sizes of 9nm, however, and larger effects may 
occur for smaller nanocrystals, or nanocrystals that are not compacted together. 
Lower-dimensional excitations were more prominent in a computational study where 
the nanocrystals were individually isolated and were very small (2-3 nm) [280]. 
In a companion experimental study, inelastic scattering from hydrogen atoms was 
measured, but the motions of hydrogen atoms in grain boundaries may not reflect the 
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Figure 28. Phonon DOS of nanocrystalline and bulk bcc 57 Fe at pressures of and 9.5 GPa, and 
hep 57 Fe at approximately 26 GPa. 

motions of the Pd atoms themselves [281]. A later molecular dynamics calculation on 
nanocrystalline microstructures of fee Ni and Cu reported similar low-dimensional 
behavior for 5 and 12 nm nanocrystals [282]. 

In a challenging experiment, Roldan Cuenya measured the 57 Fe partial DOS of 
individual nanocrystals, isolated from each other, using the NRIXS method. [283]. 
The nanocrystals of 8-10 nm diameter consisted of bcc Fe cores, surrounded by oxide 
shells, and were deposited on Si wafers. At phonon energies below 5 to lOmeV, 
exponents n in the expression g(e) oc e n were found to be 1.33, 1.84, 1.86 for three 
different samples, showing a clear low-dimensional behavior in the partial DOS of the 
Fe. This behavior depended on the shell of surface atoms of the nanocrystals, where 
the tendency toward 2D vibrational behavior was likely due to the large number of 
low-coordinated interfacial atoms. Supporting calculations showed that the lower- 
dimensional behavior extends from the surface to the core of the nanocrystals due 
to hybridization of surface and bulk modes. The evidence for lower-dimensional 
vibrations seems stronger in systems where there is a larger acoustic impedance 
mismatch between surface atoms and their surroundings, such as between isolated 
nanocrystals and vacuum. In such cases, however, it may be more appropriate to 
consider the low energy vibrations as eigenmodes of small particles [267, 268, 269]. 

In a separate NRIXS study, Roldan Cuenya, et al. studied phonons in thin, 8 nm 
layers of 57 Fe confined between layers of Cr, Co, Cu, Pd, Au, or Ag [284]. There 
was a clear suppression of high-frequency longitudinal modes of bcc Fe when it was 
confined between metals having a lower cutoff of their phonon DOS. Owing to the 
large mismatch in phonon modes, high-frequency modes of Fe cannot propagate into 
the Cu, Pd, Ag, and Au layers, and remain confined inside the thin Fe layers. Short- 
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wavelength vibrational waves in the Fe layers are expected to be suppressed near 
the Fe/M interfaces when M=Cu, Pd, Ag, and Au, consistent with observations. 
The phonon DOS of the bec Fe is shifted to lower energy, causing an increase in its 
vibrational entropy. In addition, thermal transport is also expected to be suppressed 
in these thin film multilayers. 

For nanocrystals, it is generally found that phonon spectral weight is trans- 
ferred to lower energies, causing an increase in the vibrational entropy as large as 
0.2 /ce/atom for nanocrystalline materials measured to date [273, 275, 276, 277, 285]. 
Nevertheless, it appears that the entropy of vibrations and configurations is only 
about half of what is required to give thermodynamic stabilization to a nanocrys- 
talline microstructure at reasonable finite temperatures. Further research should 
show if greater stabilization is possible by the enhancement of vibrational modes 
at low energies. Some of this increase in vibrational entropy can be offset by the 
spectral weight that appears at high energies, however. 

16.2.3. Broadening of Spectral Features of Nanocrystalline Materials 
The second prominent feature of the phonon DOS of nanocrystalline materials, 
but one that remains less well understood, is the broadening of all features of the 
DOS. This broadening is at least partially independent of the transfer of spectral 
weight to low energies. One possible explanation is that the broadening originates 
with the "damping" of the phonons, which means the rapid transfer of energy from 
one phonon state to another, likely caused by an anharmonic perturbation. Such 
behavior can originate with anharmonic interatomic potentials, as may perhaps 
exist for atoms in irregular sites in grain boundaries. These broadenings can be 
large, as shown in Fig. 27. If the effects originate with damping, the spectral shape 
of a damped harmonic oscillator may be appropriate 



where Q is the "quality factor" of the oscillator (i.e., the number of cycles for which 
the fraction e _1 of the phonons will remain) and e is the range of energies around 
the energy of the phonon, e'. The shape of D £ >(e) is similar to a Lorentzian function 
when Q is large, but D £ i(e) goes to zero as e goes to zero. For spectra measured to 
date on consolidated nanocrystals, convoluting the phonon DOS of a bulk material 
with the function D £ >(e) has proved successful in accounting for the broadening of 
the phonon DOS of nanostructured materials, except for the independent effects at 
low energies that are caused by degrees of freedom of the nanostructure itself. Figure 
29 shows such an example. 

A curious feature of the broadening of the phonon DOS of nanocrystals is that 
sometimes it can be changed without changing the size of the nanocrystallites. As 
shown in Fig. 30, the phonon DOS of nanocrystalline Ni3Fe sharpens considerably 
after the sample is exposed to deep cryogenic temperatures, even though the low 
energy part of the DOS remains essentially the same. It is possible that some of 
the atoms in anharmonic environments, such as atoms in energetically unfavorable 
sites at grain boundaries, were jostled into more favorable harmonic configurations 
by the thermal contractions induced by cryogenic exposure [286, 287, 288, 289]. 
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Figure 29. The characteristic shape of the damped harmonic oscillator function is shown in the 
center for modes at 10 and 35 meV having the quality factor Q = 7. Solid curve at bottom is the 
convolution of the damped harmonic oscillator function for Q = 7 with the experimental DOS from 
the sample annealed at 600° C (points at top). Points at bottom are experimental DOS curve from 
6nm Ni-Fe material. After [275]. 

An alternative explanation for vibrational modes at high energies was suggested 
by D. Wolf, et al. [290]. Their molecular dynamics simulations showed modes at 
energies above those of the crystallites. These high-energy modes were attributed to 
local atomic environments in grain boundaries with reduced interatomic separations. 
This explanation could be consistent with the trend of Fig. 30, assuming again that 
cryogenic exposure or low-temperature annealing were to reduce the local atomic 
distortions within grain boundaries. 



17.1. Cerium 

The T-P phase diagram of cerium metal is shown in Fig. 31. With increasing tem- 
perature, there are four solid phases at ambient pressure. Cerium is particularly 
unusual because its a and 7 phases both have the fee structure, but their specific 
volumes differ by about 17% at ambient temperature. There is a solid-state critical 
point where the two phases merge, and it is possible to take a path in T and P from 
one phase to the other with a continuous change in lattice parameter. The fcc-fcc 
transformation has attracted considerable attention, in part because of suggestions 
that the outer electron(s) of cerium become delocalized when the 7-phase collapses 
into the high-density a-phase. Nevertheless, all the phase transformations of cerium 
prove interesting, and we begin with the high-temperature phases. 

The phonon DOS above and below the j-S (fcc-bcc) phase transformation at 
999 K was the subject of two independent neutron scattering experiments [291, 292], 
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Figure 30. Phonon DOS curves at 300 K of fee NiaFe as- milled (with 6 nm grain size), before 
and after cooling to 10 K. After [275]. X-ray diffractometry showed no change in grain size after 
cryogenic exposure. 
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Figure 31. Low-pressure region of the cerium T—P phase diagram. The solid phases at zero pressure 
are: fee a, double hep /3, fee 7, bec 8. 



using two different types of instrumentation. Cerium scatters neutrons coherently, 
so it is necessary to sample a wide range of momentum transfers, Q, to get a good 
representation of the full phonon spectrum. A triple-axis experiment was performed 
with several different values of Q, chosen strategically for their sampling of the 
phonon DOS. The conclusion from this work was that the change in vibrational 
entropy was too big to explain the latent heat of the 7— # transformation. The second 
measurement with a time-of- flight instrument gave a better coverage of Q, and con- 
cluded that the change in vibrational entropy was even bigger. Both measurements 
were difficult, however, owing to the high temperature of the 7~<5 transformation, 
and the large correction required for multiphonon scattering. Giving twice the weight 
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to the second experiment, the vibrational entropy of the 7-d transformation would 
be 0.64/cB/atom. The latent heat of this transformation is 0.36A;B/atom [293, 294], 
which leaves -0.28 k-g /atom to be explained by other sources of entropy. Electronic 
entropy is the probable cause. Such a large change in electronic entropy at such a 
high temperature is surprising - one might expect the Fermi-Dirac distribution and 
other thermal broadenings of the electronic DOS to wash out any large differences 
in electronic entropy. 

The (double hcp-fcc) phase transformation is sluggish at low temperatures. 
(A 20-year study showed that 283 K is a better estimate of the transformation 
temperature than that of Fig. 31 [295].) Inelastic neutron scattering experiments 
performed on the f3 and 7 phases of cerium [296, 297] found only a small difference in 
vibrational entropy, approximately 0.09 &b / atom at the transformation temperature. 
The vibrational entropy was the largest contributor to the small latent heat at the 
transformation temperature, but other sources of entropy were much larger at lower 
temperatures. In particular, the entropy difference between the 7 and (3 phases had 
large contributions from the difference in heat capacity of the crystal field levels 
of the two phases, which was as large as -0.4 k-Q /atom at approximately 40 K. 
The difference in electronic entropy of the 7 and /3 phases made a small, negative 
contribution at all temperatures. The spin fluctuations are responsible for significant 
entropy, but curiously, the spin fluctuations in both phases are rather similar, and 
do not contribute a significant entropy to the phase transformation. 

The celebrated a-7 (fcc-fcc) phase transformation in pure Ce is usually induced 
by pressure, but it is possible to obtain this phase transformation as a function 
of temperature in the alloy Ceo.gTho.i. An inelastic neutron scattering experiment 
reported differences in the phonon DOS of the a- and 7-phases, but these were 
responsible for only a negligible difference in the vibrational entropy of the two 
phases [298]. The entropy of the a-7 transformation in Ceo.gTho.i originates almost 
entirely from the entropies of electrons and spins. These contributions are shown in 
Fig. 32. The electronic entropy from disorder in the electron occupancies near the 
Fermi surface makes a modest contribution to the difference in entropies of the two 
phases. The ground state spin fluctuation entropy was deduced from the quasielastic 
magnetic scattering, and the crystal field entropy from the inelastic excitation of a 
broad crystal field transition at 16 meV. These magnetic contributions are very 
large, and are the main source of the entropy difference between the two fee phases 
of Ce . 9 Tho.i. 

The a-7 transformation in pure Ce seems to have a large change in vibrational 
entropy, unlike the a-7 transformation in Ceo.gTho.i. Jeong, et al., performed diffrac- 
tion measurements on pure Ce metal, crossing the a-7 phase boundary by varying 
T at fixed P, and by varying P at fixed T [300]. They observed a discontinuity in 
the Debye- Waller factor at the phase transformation, indicating that the vibrational 
entropy of 7-Ce is 0.75±0.15 &B/atom larger than that of a-Ce. The mean-squared 
displacements of atoms in a phonon of frequency co go as w~ 2 , so the Debye-Waller 
factor emphasizes low- frequency phonons more strongly than the — In uj for the 
vibrational entropy of Eq. 16. This difference could be a problem if the shape of 
the phonon spectrum differs between different phases, as is the case for Ceo.gTho.i. 
Nevertheless, Jeong, et al., performed diffraction measurements on this alloy, too, 
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Figure 32. Contributions to the entropy difference between the 7 and a phase of Ceo.9Tho.1- The 
thick curves is the sum of the four components. The electronic component was adjusted so that the 
sum equaled the entropy measured by [299] . 

and found small changes in the Debye- Waller factor across the phase transformation 
[301], qualitatively consistent with prior neutron scattering results. It seems that 
pure Ce has a large vibrational entropy of its a-7 transformation (accounting for 
about half the transformation entropy), but the addition of 10% Th suppresses the 
vibrational entropy of the a-7 transformation. 

17.2. Uranium and Plutonium 

Somewhat like cerium, uranium assumes 3 solid phases at different temperatures 
and ambient pressure. At low temperatures, the orthorhombic a-phase shows a 
number of unusual properties such as a charge-density wave transition [302, 303]. 
Such an instability in electronic structure suggests that the atom displacements in 
phonons may alter the electronic structure, and the interatomic forces and phonon 
frequencies may change with temperature. A temperature-dependent electronic ef- 
fect on phonons was in fact found in a neutron inelastic scattering study on a-U 
[304]. The quasiharmonic phonon softening in a-uranium expected with Eq. 83 
and the approach of Section 6.3 is small, owing to the low bulk modulus. The 
actual phonon softening to a temperature of 940 K gives an increase of vibrational 
entropy of 1.5/cB/atom, approximately a factor of 10 larger than the quasiharmonic 
prediction. Typical anharmonic potentials cannot account for such a large softening 
of the phonons. The clear implication is that the interatomic forces are weakening 
with temperature. It is not known if this originates with thermal excitations of 
electron energy levels, electron energy level broadenings caused by shortened electron 
lifetimes, or by changes in electronic levels with temperature that are caused by 
larger atomic mean-squared displacements, for example. 

At ambient pressure, plutonium has even more solid phases, six, at different tem- 
peratures, having atomic volumes differing by as much as 20%. With such differences 
in local atomic distances, the 5/ bands of Pu differ considerably in their widths for 
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the different phases. For the denser phases, the Peierls distortions of the 5/ bands are 
responsible for the formation of crystal structures of low symmetry [305, 306]. The 
fee (5-phase of Pu is intermediate in its 5/ localization, and temperature causes some 
odd behavior, such as negative thermal expansion. The fee <5-phase can be stabilized 
over a broad range of temperature by alloying with 2 at.% Ga or 5 at.% Al (although 
the thermal expansion becomes positive for the latter). Thermal neutron scattering 
from Pu is possible for samples with a high enrichment of the Pu isotope, which 
has an acceptable absorption cross section. It is also possible to measure phonon 
dispersions with inelastic x-ray scattering, which does not require isotopic control, 
and a full set of phonon dispersions were measured on a small single crystal of Pu- 
2 at.% Ga at ambient temperature [307]. An anomaly in the low transverse [Oil] 
phonon branch, and small anisotropic elastic constants were attributed to electron- 
phonon interactions. Measurements of the phonon DOS at temperatures of 300 K 
and below were performed on 35 g of 242 Pu-5at.% Al [308]. There was generally good 
agreement with the results of [307], but from 27 to 300 K there was a large overall 
softening and broadening of the phonon spectrum. Similar to the case of U, for 
Pu the quasiharmonic softening expected from Section 6.3 was more than a factor 
of ten too small to account for the phonon softening - the softening occurs with 
relatively little change in lattice parameter. Such behavior seems plausible because 
temperature should have large effects on electrons in very narrow 5f bands, but 
detailed predictions of electron-phonon interactions in Pu are only now becoming 
available [309]. 

The low temperature phases of Pu tend to have smaller specific volumes than 
the high temperature phases, and the actinide metals to the right of Pu (Am, Cm, 
Bk...) tend to be much larger in specific volume. It has been suggested that for 
the larger phases of Pu, and for the heavier actinides, the 5/ electrons become 
localized and non-bonding [310]. With these changes we may expect lower phonon 
frequencies. Wallace used calorimetric data to show that this trend is generally 
true, and also showed that vibrational entropy plays a big role in stabilizing the 
high temperature phases [311]. The low temperature a-phase of Pu was also shown 
to be rather anharmonic. 
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III. Vibrational Entropy of Materials at Higher 

Temperatures 



18. Interactions between Thermal Excitations of Electrons and 

Phonons 

18.1. Overview of Thermal Effects on Electron and Phonon Energies 

The full Hamiltonian of the crystal, H, is 

H = H n + H c + H cp (140) 

where H n is the contribution from the nuclei (i.e., phonons as in Eqs. 23 or 141), H e 
is the contribution from the electrons (including electron-electron interactions and 
thermal electronic excitations), and H ep is the contribution from electron-phonon 
interactions (EPI). Figure 33 attempts to illustrate some concepts and definitions 
of high temperature phonon thermodynamics, especially the adiabatic EPI. Six dia- 
grams are shown, with low temperature effects on the bottom and high temperature 
effects at the top. Phonon concepts for a fixed interatomic potential 5>(r) are on the 
left; electron concepts with a simple free electron band structure are on the right. 

At low temperatures, the vibrating atoms explore a small part of their interatomic 
potential wells. A parabolic fit to the interatomic potential is appropriate (at least 
without zero-point effects), and the harmonic model follows. At higher temperatures, 
the "quasiharmonic" approximation is often used, where Fig. 33 shows that the 
parabola is altered in position and curvature to accommodate small shifts in the 
effective potential as the atoms vibrate about their points of equilibrium. (The 
thermal spread in potential energy, of order k^T, is shown by light shading.) The 
quasiharmonic model assumes that phonons are eigenstates of a phonon Hamiltonian 
for the expanded volume of the crystal, and it does not predict any shortening 
of phonon lifetimes. This is less appropriate at higher temperatures, where the 
cubic and quartic contributions to the potential energy are needed to understand 
the reduced phonon lifetimes, and shifts with temperature that may not depend 
quasiharmonically on the volume of the crystal and its bulk modulus (Section 6.3). 
"Anharmonic" corrections to the quasiharmonic model are needed to account for 
"phonon-phonon interactions" (PPI). 

The electron states at right in Fig. 33 are filled to the Fermi level at T = 0, 
but at higher temperatures the shading illustrates a thermal spread of energies 
associated with the Fermi-Dirac distribution. At higher temperatures, as the nuclei 
move further about their equilibrium positions, the electron band structure is al- 
tered. If this occurs without thermal excitations of electrons, the change in energy 
is that for anharmonic phonons, already included in H n . On the other hand, if 
there are atom displacements from many phonon excitations plus there are thermal 
electron excitations, there is a change in energy that depends on both excitations. 
The conditions for obtaining this extra energy are illustrated at the top of Fig. 33. 
This is an adiabatic electron-phonon interaction, and makes up H cp of Eq. 140. 
Although differing in its thermal factor, the adiabatic EPI has similarities to the 
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Figure 33. Effects of temperature, increasing from bottom to top, on phonons (H n , left), electrons 
(H c , right), and the existence of the adiabatic EPI when phonons alter the electronic bandstruc- 
ture in the presence of electron excitations (H Bp , top). Increasing width of light shading indicates 
increasing thermal energy spreads. 



nonadiabatic EPI responsible for Cooper pairing in conventional superconductors, 
and until recently it was assumed to be a weak phenomenon that could be neglected 
at elevated temperatures. 

18.2. Formalism of Phonon-Phonon Interactions 

To understand phonon-phonon interactions, we consider in detail the Hamiltonian 
for nuclear motions H Q , extending the Taylor series of Eq. 23 



2m 2! 

na K ' a i 



+ oT X] ^oca'a" kk'k" U aK U a ' K > U a " k" 



Ka K ' a i K » a " 



+ JT^X X X! ®aa'a"a"'K.K<K"K"< U aK U a / K > U a " K " U a ">K"> + ••• j (141) 



4! 

™ re' a' K"a" K"'a"' 



where we have simplified the development to one atom per unit cell, so u aK is the 
displacement of atom k along Cartesian direction a. 

The transformation to normal coordinates for an infinite periodic crystal is 



Ur = -L= Y u r - , (142) 



and by Fourier inversion 



Ur = -i= T Ur e +i ^ . (143) 
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Substituting Eq. 143 into Eq. 141, we identify the Fourier transform of $ aa ' KK ', 
and use it to define the dynamical matrix of Eqs. 30 and 31. Our interest now 
is in the higher-order terms of Eq. 141, which cause anharmonic behavior when 
the atoms are vibrating with substantial amplitude. Using Eqs. 141 - 143 we can 
assess phonon-phonon interactions with a classical approach, which is valid at high 
temperatures. 

There is, however, a clarity to be gained in working with the quantized phonon 
field (a displacement field in the crystal) , especially when considering excitations of 
nonlocal electron states simultaneously with phonon excitations. The phonon field 
operator is 

H=% +a U= A U- (144) 

The raising and lowering operators at and ar create and annihilate the phonon fcj. 

hi i 

They are constructed from the momentum and position operators of the Hamiltonian 
when a plane wave solution is considered for one ki. The reader should have refer- 
ences showing the properties of and a for a simple harmonic oscillator [312], and 
the at and ag work the same way for a single normal mode of a harmonic crystal. 
The formalism is now in the "occupation number representation," where n phonons 
of wavevector ki are created with n raising operations as (o^j |0)g = (n!) _1 / 2 |n)g. 
The phonon field operators are related to the displacements (now operators) as 
Fourier coefficients 



u 



^=?V^7* )eiE,F "^- (145) 



Note the similar roles of the phonon field operator A? in Eq. 145 and the vector 

Ki 

Ur in Eq. 143. 

Ki 

By substituting Eq. 145 into Eq. 141, the dynamical matrix is recovered from 
the quadratic term, with prefactors such as h/(4u>^ u>^ ). We define similarly the 
transformed higher-order terms, so Eq. 141 becomes 

hi 

+ E E E kj,k k ) A %i A h A %k 

hi h j h 

+ E E E E V&, kj,k k , h) A k A %j A %k A %i + ... , (146) 

hi hj hfc h^ 

where the first sum on the right-hand side includes both the kinetic energy and 
the harmonic part of the potential energy (following standard development [313]), 
and the phonon number operator is a|- ar . The V are related to the $ as Fourier 
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transforms. For example, 



-+ - - 1 1 ffr\ 3 1 - 
V(ki,kj,k k ) = — — — 3 - 5(ki + kj + k k - g) 



X 



]T e Q (^)e a ,(fe j )e Q »(4)e +i(fcl+ ^ +fcfc) ^ , (147) 

{<*}{?*} 



where the factor 5 signifies that the sum k{ + kj + k k equals a reciprocal lattice 
vector, g. For terminology, the terms with V(ki, kj, and V(ki,kj,k k ,ki) in Eq. 
146 give the energies of "phonon-phonon interactions" because they alter the phonon 
energies when larger vibrational displacements are present in the crystal. 

Details of the phonon dispersions are immediately relevant when attempting to 
calculate anharmonic behavior, and several concepts have been developed in detail. 
For example, if the phonon dispersion relation e(k) were linear, it would be easy 
to add two wavevectors k + k' = k" , and simultaneously conserve energy e + e' = 
e". The number of three-phonon processes increases as k 2 , at least for small k. 
Many phonon dispersions e(k) are concave downwards, however, so simultaneous 
energy and momentum conservation is often impossible unless at least one phonon 
lies on a different branch. Phonon processes that satisfy momentum and energy 
conservation depend on the symmetry of the dispersion relations, and on the crystal 
structure. For larger phonon wavevectors k, momentum conservation is possible by 
adding a reciprocal lattice vector. The idea is that the momentum transferred to the 
entire crystal occurs with zero energy because of the large mass of the crystal. Such 
"umklapp" processes allow many more three-phonon interactions, but the phonon 
wavevectors must be of length comparable to the reciprocal lattice vector if these 
vector additions are possible. 

An important question is how many terms are needed in Eq. 146 to account 
accurately for anharmonic behavior. An approximate theory of the lifetimes of 
optical phonons was developed with 3-phonon processes [314], and was extended to 4- 
phonon processes to include energy shifts [315]. This approach has found application 
for the analysis of the thermal shifts and widths of Raman peaks, sometimes with 
applications to thermal conductivity (e.g. [316]). We expect the higher-order terms in 
Eq. 141 to be progressively smaller, but the requirement of energy conservation and 
the ^-function for momentum conservation in Eq. 147 restricts the allowable three- 
phonon processes, owing to restrictions imposed by the phonon dispersions e(k). 
Fourth-order anharmonicity can therefore play a role as significant as the third. 
Another complexity is that four-phonon processes can be generated from three- 
phonon processes. The idea is that two three-phonon process such as k + k! — ► k" 
and k" — > k!" + k"" can be combined through intermediate states {/?'} to become 
a four-phonon process in second-order perturbation theory, which involves matrix 
elements as E^'" ' ,k""\H"'\k"){k"\H"'\k,k'). With effort, this field theoretic ap- 
proach has been used for some physical predictions as described below. Approximate 
summations of the perturbation series have been performed, and some results for 
anharmonic heat capacity were reported for solids with Lennard-Jones potentials 
[318]. This approach also has conceptual value for understanding the meaning of 
anharmonicity. 
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18.3. Formalism of Electron-Phonon Interactions 

The large difference in energy scales of phonons and electrons can motivate the 
separation of the Hamiltonian of the solid into a term with nuclear coordinates for 
the phonons H n , and a term with electron coordinates H e (Eq. 140). The energy of 
the crystal deformation caused by a phonon originates with the electrons, of course, 
but although this potential energy of deformation is electronic in origin, it transfers 
to kinetic energy in the motion of the nuclei. 45 If the electrons were always in their 
ground states, all this energy of deformation would be associated with H n alone 
(i.e., the phonons). Treating the electron system and the phonon system as two 
independent thermodynamic systems becomes inconsistent at finite temperature, 
however, because the presence of phonons alters the thermal excitations of electrons. 

Very generally, the electron-phonon interaction (EPI) requires the coordinates of 
the electrons {r^ 1 }, and coordinates of the nuclei {r^} 

tfep = £V^lV7). (148) 

First assume the nuclear motions are slow enough so the electron levels adapt 
continuously to the evolving structure, i.e., the Born-Oppenheimer approximation. 
This "adiabatic approximation" does not allow the nuclear kinetic energy to alter 
the electron states. 46 The adiabatic EPI increases with the number of phonons, 
but it also scales with the thermal excitations of electrons. (Effects from how the 
interatomic potential depends on phonon populations alone are already counted as 
phonon-phonon interactions.) 

In this adiabatic picture, unique electron eigenstates exist for a snapshot of the 
nuclear thermal displacements, but these evolve continuously as the nuclei move. 
Taking for reference a ground state with the nuclear positions of a perfectly periodic 
crystal and the electronic structure at K, we seek to account for how the displaced 
nuclear positions alter the thermally-excited electronic energy levels (changing the 
electronic excitations induced by temperature). For a given electronic structure, the 
adiabatic electron-phonon interaction is proportional to the number of phonons, 
n(e, T) + 1/2, times the difference of electron occupancy with respect to the ground 
state, f(T) — /(0), where f(T) is the Fermi-Dirac distribution. When there are 
sharp features in the electronic DOS near the Fermi level, there will be a shift of 
electron chemical potential with temperature, and other temperature-dependences 
associated with the thermal sampling of the excited electron states. 

How the electron eigenstates change with nuclear displacements is a problem that 
generally requires considerable detail about the states of electrons and phonons. A 
simple approach is to consider an electronic band formed from isotropic s-electrons, 
and a uniform dilation as may be associated with longitudinal phonons of long 

45 For a harmonic phonon, the electron-phonon interaction affects equally the energies of the 
electrons and the phonons. 

46 The "non-adiabatic" electron-phonon interaction accounts for how the electronic states are 
altered by the nuclear kinetic energy, not the potential energy of displaced nuclei as for the adia- 
batic case. The non-adiabatic EPI does not require thermal activation, and can be responsible for 
superconductivity. 
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wavelength, 

H»=-J2DA(rt). (149) 

A 

where A(r| 1 ) is the fractional change in volume at rf, and D is a "deformation 
potential," typically a few eV. 47 This simple approach is convenient because the 
A(r^) is related to the displacement u as through its divergence: A(rf) = V ■ u, so 

tf e D p = -££V-u(rt). (150) 

A 

For longitudinal phonons 

e{h) = £ , (151) 

and by taking the divergence of Eq. 145 we have, now in the occupation number 
representation 



Hg = -iDyyJ—^ \h\ e^ 1 A? (152) 

rf ki V *i 

Consistent with the second quantization formalism for phonons, fermion field 
operators are used for the electrons 



which when applied to a multi-electron state, places an electron in the state k\ (the 
operator removes an electron from this state). When small, in Eq. 152 

is a perturbation that mixes electronic states as f\H^ p \^i). This integral can be 
evaluated by using the fermion field operators C and C\ where the exponential phase 
factors cancel the /c-space integration unless there is a conservation of wavevector. 
This forces the same total k\ + ki for the creation operators as for the annihilation 
operators 



/ Z electron phonon j , 

Hg = -wJ-£- y y J^-AgCt -C E , (154) 
plus an analogous term with Al CV , r Ci for phonon creation. For Eq. 154 we also 

ki k\ 

used a linear dispersion relationship coLih) = cjjzi, appropriate for long- wavelength, 
longitudinal acoustic phonons. In accounting for the effects of phonons on the elec- 
tron states, it is usually necessary to consider a more general case than the isotropic 
deformation potential of Eq. 149 as 



I oc (k x + k x 



e^)-V^(rt, Tj) &a) j (155) 



47 It is not typical to have electron energies shifted by several eV under high compressions, so the 
D typically pertain to small elastic distortions over a finite volume of crystal. 
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where the gradient is with respect to the nuclear coordinates. 8 

In general, there are two lower order terms that are used to describe the electron- 
phonon interaction. The first is a generalization of the previous result Eq. 154 

4p = EE ycv ^ HF\ X ^ ■ (we) 

It accounts for the processes where an electron is excited from state k\ to k\ + 
ki, simultaneously with the annihilation of a phonon in state ki, or creation of a 
phonon in state —ki. Momentum conservation restricts significantly these three- 
phonon processes. The second low-order term for electron-phonon coupling is H 2 p , 

and it includes a factor from second-order perturbation theory, J2p(k\^\k r )(k f \fi\k). 

Note that the initial and final states of the electron have the same k, so there is no 
change of electron state by this process. Nevertheless, the scattering into the virtual 
states {k 1 } serve to alter the self-energy of the electron. 

For calculating the thermodynamic effects of electron-phonon coupling, it is more 
useful to know how the electron-phonon coupling depends on energy, rather than 
momentum as in Eq. 156 for . This requires averaging over all electron states 

near the Fermi surface separated by kj and differing in energy by a selected fku. 
Most work on electron-phonon coupling has focused on superconductivity, where 
hio is rather small, and the electron states are close enough to the Fermi surface 
that it is reasonable to use ground state Fermi surface properties. The Eliashberg 
coupling function a 2 g(uj), where g(u>) is the phonon DOS, accounts for all scattering 
between electron states on the Fermi surface. The function a 2 g(u>) often looks similar 
g(uj), although some parts of the spectrum are given different weights. An important 
moment is calculated by weighting by and can be shown to be 

Jo uj M{uj 2 ) v ' 

where (I 2 ) is an average of 1*1 of Eq. 155 (or its analog) over all directions, (lo 2 ) 
is a second moment of the phonon DOS, N(Ep) is the electron DOS at the Fermi 
surface, and A is the "electron-phonon coupling parameter." The strength of the 
electron-phonon coupling scales with A. 



19. Phonon-Phonon Anharmonicity 

19.1. Damped Harmonic Oscillator 

The equation of motion for a classical damped harmonic oscillator is 

d^x dx 

m— + b— - kx = . (158) 



"Frolich coupling" is another example of an EPI, where longitudinal optical phonons create an 
electric field in an ionic solid. In this case, the Fourier transform of a Coulomb potential gives a 
factor of Ank^ 2 , and with two ions involved, goes as fcr 1 (instead of the for H® p in Eq. 

154). Transition rates, or scattering rates, go as the square of these quantities. 
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We define the dimensionless quality factor Q as 

a s ^ . 

A solution of Eq. 158 is 

x(t)=X e-^ 2Q ^ 1 ^ t , (160) 

where the frequency u> differs somewhat from that of the undamped harmonic 
oscillator (for which uq = ^Jk/m) 



U = UJ «i l + (2Q)-* - (161) 



When Q is large, the frequency shift A is approximately 

UJQ ( 1 \ 2 



■ (162) 

The energy spectrum of the damped harmonic oscillator was given in Eq. 139 
for e' = fouj. When Q is large, Eq. 139 has the form of a Lorentzian function with 
half-width-at-half-maximum, T, of 

r-£. (163) 

For a damped classical harmonic oscillator, the ratio of peak broadening to peak 
shift is 

J = -4Q . (164) 

Equation 164 is convenient for testing against experimental measurements. 

When the vibrational entropy is evaluated with Eq. 14, there is a change to S v fo 
when a DOS curve g(e) is broadened by the damped harmonic oscillator power 
spectrum of Eq. 139. The shape of the damped oscillator spectrum is asymmetrical, 
going to zero at e = 0, and having more weight on its higher energy side. At high 
temperatures, where the phonon occupancy is large, the broadening serves to reduce 
the vibrational entropy. 49 (There is a compensating effect from the lowering of the 
frequencies by damping, but this effect is smaller.) 

The damping of a classical harmonic oscillator is understood as a dissipative 
process that transforms the energy of a macroscopic system into heat. Phonons arc 
the heat, however, so the damping of a phonon cannot be understood in the same 
way. For example, adding a dissipative term to the Schrodinger equation would 
dissipate probability, not energy. It is necessary to treat phonon lifetimes as the 
transfer of energy between phonon modes, so a multi-body quantum mechanics is 
immediately necessary. The coupling of anharmonic quantum oscillators allows prob- 
ability to transfer between the different oscillators, sometimes with energy shared 



49 When two DOS curves such as g a (e) and g^{e) of Eq. 16 are both broadened by damped 
oscillator functions of the same Q, there is no effect on AS^ a . 
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by more than one. Such a model is convenient because it conserves energy, and a 
Hamiltonian can be written for the problem. If only a few oscillators are involved, 
however, a state will recur with some characteristic time. This may or may not be an 
issue for understanding experimental measurements of phonon lifetimes, but it will 
lead to confusion about the meaning of dissipation and irreversibility. The solution 
to the problem is to have a dense distribution of oscillators in the heat bath, which 
is an appropriate model for a macroscopic system with many degrees of freedom. 
For this reason, the damping of a quantum oscillator is sometimes analyzed as the 
imaginary part of the dynamical susceptibility of an external field interacting with 
the oscillator [319]. 

Self-consistent phonon theories offer an alternative to the perturbation theory 
approach described in Section 19.2. In essence, these methods make a modification 
to the dynamical matrix of harmonic theory. Horton and Cowley develop a simple 
one-dimensional example where the force constant, d 2 V/dx 2 , is convoluted with a 
Gaussian function in x [320]. A set of equations is solved for consistency between 
the breadth of this Gaussian and the phonon frequency. This first-order approach 
uses a symmetrical interatomic potential without odd terms in x. This is also true 
for the more general first-order self-consistent phonon theory developed by Wallace, 
who then demonstrates interactions between the first-order phonons from the cubic 
terms in the potential [46]. 

19.2. Perturbation Theory 

The usual approach to calculating how the frequencies and lifetimes of phonons 
change with temperature is with perturbation theory, using harmonic phonons as a 
basis set, and the last two terms of Eq. 146 as perturbations that alter the energies of 
the vibrational states. The phonon lifetime is calculated from Fermi's second golden 
rule, where transition rates are calculated individually by expanding the operators 
A? in Eq. 146, and using these to couple the initial and final states of harmonic or 
quasiharmonic phonons. 

Perturbation theory is performed for the atoms about their equilibrium config- 
uration at elevated temperature, i.e., the thermal expansion consistent with quasi- 
harmonic theory is already taken into account. For centro-symmetric crystals, the 
cubic terms in Eq. 141 with ^ a a.' a" kk 1 k" u an « 'k' u o"k" vanish because they are odd 
in the atom displacements, u. For example, the cubic energy contribution from the 
positive displacement of an atom with respect to its neighbor on the right is exactly 
cancelled by its negative displacement with respect to the atom at its left. For 
quantized vibrations, there is no energy change when V(ki, kj, kk) from Eq. 146 (or 
147) is treated as a first-order perturbation of quasiharmonic phonon states. The 
first non-vanishing first-order perturbation is the fourth-order term V(ki,kj,kk,ki) 
in Eq. 146. It is physically reasonable, but argued more rigorously by Van Hove [321], 
that with respect to the harmonic term of relative strength 1, the energy of the n th - 
order perturbation goes as (u/a) n ~ 2 , where u is the displacement, and a is the INN 
distance. We note that the fourth-order term in first-order perturbation theory is 
of relative strength (u/a) 2 , but the third-order term in second-order perturbation 
theory is of equal strength, i.e., (u/a) 3 " 2 x (u/a) 3 " 2 = (u/a) 2 . The standard treat- 
ment of anharmonicity by perturbation theory includes these two parts, the quartic 
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Figure 34- Summary of the eight terms from the cubic term of Eq. 146 obtained by expanding the 
operators = a%. + a 1 g • 

term in first-order and the cubic term in second order. Many individual terms must 
be calculated, however, as summarized in Fig. 34 for the cubic term of Eq. 146. Not 
all of these terms conserve energy, such as the top and bottom processes in Fig. 34. 
They nevertheless do need to be considered as intermediate states (virtual states) 
in the second-order perturbation theory treatment of the cubic term. 

Maradudin and Fein [322] obtained expressions for the phonon frequency shift in 
the presence of third- and fourth-order anharmonicities of Eq. 146, and these results 
have been discussed in a number of other references [99, 323, 324]. The phonon 
frequency shifts A(kj,ujj,T) and damping T(kj,u)j,T) have the following features 

— The quasiharmonic contribution to A(T) is identified. (It is zero if the volume 
of the crystal is fixed.) 

— For the phonon lifetime broadening owing to scattering between quasihar- 
monic phonon states, the most important term is from the cubic anharmonicity. 
This term appears in first order perturbation theory, so it includes only those 
processes in Fig. 34 for which energy is conserved. 

— For the shifts in phonon frequencies, the lowest-order anharmonic terms are 
the cubic term in second-order perturbation theory, 50 and the quartic term 
in first-order perturbation theory. These can be comparable in strength, and 
comparable to the quasiharmonic part. It is not obvious what sign is expected 
of these anharmonic terms. They are nonzero even when there is no thermal 
expansion. 

— At high temperatures, A(T) and T(T) have contributions proportional to both 
T and T 2 . These usually do not reconcile with the classical prediction of Eq. 
164, in part because Q typically increases with temperature. 

It is expected that the results for A(T) and T(T) from perturbation theory are 
reliable for small atom displacements, u/a. The practical range of temperatures 



In the special case when three times the phonon wavevector equals a reciprocal lattice vector, 
the cubic term can alter the phonon energy in first order. 
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where these results are applicable is still a topic for research, but it is known that 
these results are not reliable near the melting temperature of rare-gas solids [325], 
and are not expected to work for quantum solids such as helium at any temperature. 
Unfortunately, the effort to extend these results to higher order seems very difficult, 
although there has been an impressive success in summing diagrams for a well- 
behaved model of a rare-gas solid [326]. 

19.3. Anharmonicities of Individual Phonons 

Experimental measurements of the temperature-dependence of phonon frequencies 
and lifetime broadenings have been used to assess the relative importances of the 
cubic and quartic terms of the interatomic potential. These measurements have been 
performed on insulators such as alkali halides, for which there is no role for thermal 
excitations of electrons or electron-phonon interactions. A number of studies have 
used optical spectroscopies, which sample phonons of small k. For example, for the 
transverse optical mode at k = 0, the cubic and quartic anharmonicities for AgBr 
are negative, but positive for TIBr [327]. Neutron scattering studies of phonons at 
the Brillouin zone boundary in AgGaSe2 with the chalcopyrite structure showed 
that the anharmonic contributions dominated for an optical phonon, whereas the 
quasiharmonic approximation was adequate for an acoustic phonon [328]. 

An anomalously-large phonon damping was found at high temperatures in CaF2, 
a superionic conductor [329]. The phonon frequencies were approximately consistent 
with effects of volume expansion, but the Griineisen parameters for different phonons 
varied from less than 1 to greater than 3, and showed some temperature dependence. 
Especially large non-harmonic effects were found for a mode involving the F~ ions. 
It might be argued that CaF2 is approximately quasiharmonic, but the Griineisen 
parameters depend on T, and the large dampings of modes are inconsistent with the 
quasiharmonic model. A recent Raman study of Hf02 showed a similarly large range 
of mode Griineisen parameters, and phonon lifetime broadenings that were approx- 
imately linear with temperature [317]. By comparing to computational results, the 
largest anharmonic effects could be attributed to modes dominated by the motions 
of oxygen ions. 

Phonon DOS curves were obtained at several temperatures for both the hP24 
and fee phases of C03V [121], giving information about how the different phonon fre- 
quencies shift with temperature. An average shift can be used to obtain a Griineisen 
parameter 7, and with a thermal expansion measured by diffractometry, a quasi- 
harmonic prediction for the difference in phonon DOS between temperatures T and 
T + AT is 



For the two phases of C03V, Eq. 165 gives the two dashed curves in Fig. 35. The 
solid curves are the actual differences in phonon DOS curves measured at two tem- 
peratures. For the hP24 phase, the agreement is perhaps qualitatively correct, but 
for the fee phase the agreement is poor. Some ad-hoc adjustment of the Griineisen 
parameter for the fee phase is helpful but not sufficient. For both phases of C03V, 
the low transverse acoustic modes are more temperature dependent than the others, 




(165) 
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Figure 35. Top: Phonon DOS curves of the fee and hP24 phases at high temperatures. Bottom: Solid 
curves labeled AhP24 cx and Afcc cx are differences between pairs of phonon DOS curves measured 
at two elevated temperatures (for hP24 at 1308 and 1073 K), (for fee at 1513 and 1333 K). Dashed 
curves labeled AhP24 ca i c and Afcc ca i c were calculated by assuming the same constant Griineisen 
parameter for all phonons with Eq. 165. After [121]. 
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Figure 36. Phonon frequency shifts calculated with a Born-Mayer potential fit to give the correct 
phonon frequencies of Zr at high temperature. After [73]. 



indicating that they are beyond the predictions of quasiharmonic theory, whereas 
the longitudinal modes seem more quasiharmonic. 

Figure 36 shows anharmonic phonon frequency shifts for bec Zr calculated with 
perturbation theory for phonons along different crystallographic directions. Approx- 
imately, the amount of shift scales with the energy of the phonon, but this varies for 
the different phonon polarizations and is not an exact relationship. Notice, however, 
that towards the V point (i.e., k = 0), the shift goes to zero, as does the phonon 
energy. There are some discrepancies between the results of Fig. 36 by May, Miiller 
and Strauch [73], and earlier theoretical work by Ye, et al., [103] on bec Zr. 

Ye, Ho, Chen, and Harmon evaluated phonon-phonon anharmonicity in Mo by 
first-principles frozen phonon calculations. The energies from their density functional 
theory calculations were used to estimate central force and bond-bending contribu- 
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tions to the third and fourth-order anharmonicities of Eq. 146 [104]. The phonon 
frequency shifts and phonon broadenings were then calculated by the perturbation 
theory discussed in Section 19.2. An interesting finding was that a central-force 
model of the anharmonic contributions was inadequate for bec Mo, although it was 
successful in their earlier work on bec Zr [103]. The authors noted that the central 
force anharmonicity tended to cause phonon stiffening with increased temperature, 
whereas the bond-bending anharmonicity tended to cause a softening. Even with 
the addition of bond-bending forces, however, it was found that there was a dis- 
crepancy at the (100) H-point when the results were compared to the experimental 
measurements of Zaretsky, et al., [102]. This direction corresponds to a Fermi surface 
nesting, and a Kohn anomaly (discussed in Section 20.1.3). 

Phonon linewidths are dominated by 3-phonon interactions in first order, specif- 
ically the two middle processes shown in Fig. 34. Some special combinations of 
phonons interact strongly and decay quickly. For example, in gallium phosphide 
(GaP) the transverse optical (TO) phonon with k = can decay into a longitudinal 
acoustic (LA) plus a transverse acoustic (TA) phonon with opposite momenta, and 
the energy is matched to that of the original TO phonon. For GaP under pressure, 
the phonon lifetime broadening was probed by Raman spectrometry [330, 331] The 
role of pressure was to suppress the matching of energy and momentum for this 
TO^LA+TA process, narrowing the Raman linewidth of the TO mode. The pos- 
sibilities for matching polarizations, energies, and momenta in 3-phonon processes 
can be complicated, with special combinations giving "Fermi resonances" such as 
this one. 

Calculating the /c-dependence of phonon lifetime broadening generally requires 
rather accurate information about the individual phonon energies. Simple rules 
for anharmonic phonon lifetime broadening are not yet evident, nor are rules for 
anharmonic phonon frequency shifts. 

19.4. Anharmonicity and Heat Capacity 

Another approach to finding anharmonicity is by assessing heat capacity data. 
Although heat capacity can include contributions from electrons and spins, for 
nonmagnetic insulators such as oxides at high temperature, it is appropriate to assign 
all the non-harmonicity beyond quasiharmonic behavior to phonon anharmonicity. 
This was the approach used by Anderson [332], who compiled an appendix on high 
temperature thermophysical properties of a number of minerals. Two minerals show- 
ing substantial anharmonicity are Mg2Si04 and Grossular garnet (CaaAk^SiO^)). 
Pyrope garnet and an olivine mineral showed modest anharmonicity, whereas AI2O3, 
MgO, KC1, CaO did not show much deviation from expected quasiharmonic behavior 
[332]. Figure 37 shows some of these reduced heat capacity data at high tempera- 
tures. The curve for Grossularite shows a particularly striking increase above the 
harmonic curve (and the classical limit), but it may indicate a change in defect 
density or another structural change in the material. Studies of phonons in FesSi 
alloys at temperatures of 20, 930, and 1100°C showed large phonon broadenings, 
approximately proportional to phonon wavevector [333]. These broadenings were 
interpreted as effects of chemical disordering in this material, which undergoes an 
order-disorder transformation at elevated temperatures. 
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Figure 37. Heat capacity at constant volume Cv for three minerals, and a hypothetical Debye solid 
with a Debye temperature of 700 K. Data reduced from tables in Ref. [332]. 

Wallace and colleagues [50] performed a similar thermodynamic analysis for 
metals, and isolated anharmonic contributions to the entropies. They obtained the 
quasiharmonic contribution from the effects of thermal expansion on the phonon en- 
tropy and on the electronic entropy. Eriksson, Wills, and Wallace performed density 
functional theory (DFT) calculations of electronic structure that allowed them to 
account for the electronic entropy by calculating the electronic DOS at T = 0, and 
using the Fermi-Dirac distribution to sample the electronic excitations across the 
Fermi surface [334]. After subtracting the quasiharmonic and electronic entropies 
from the total calorimetric entropies obtained from the literature, the remainder 
was the anharmonic entropy attributed to phonon-phonon interactions. This S a (T) 
is shown in Fig. 38. Differences in S a (T) for fee and bee metals are evident in Fig. 38, 
and supported by further results tabulated in Ref. [50]. The anharmonic entropies of 
fee metals tend to be similar, being relatively small and negative. For elemental bec 
metals, the anharmonic entropy can be large. There is significant variation in the 
anharmonic entropies of different bee metals - Fig. 38 shows substantial differences 
between Mo and Ta, for example. Heat capacities and other thermodynamic data 
of Mo and W were assessed by Fernandez Guillermet and Grimvall, who reported 
similarly-large anharmonicities in these bee metals [335, 336]. 

Vibrational entropy is fundamentally determined from the volume of phase space 
in position and momentum coordinates that is explored by the atoms in the ma- 
terial. The phase space explored by coupled anharmonic oscillators is much less 
well understood than for harmonic solids, and has produced modern research into 
nonlinear phenomena such as chaos. Nevertheless, it is hoped that the wealth of 
results on harmonic solids can be adapted to the anharmonic solids. One approach 
is to use the energies of phonon creations and annihilations measured by inelastic 
scattering experiments - although the measured phonon spectra show shifts and 
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Figure 38. Anharmonic entropy for bec and fee transition metals, scaled to the melting temperature 
T m . After [334]. 

lifetime broadenings, the experimental phonon DOS measured at high temperature 
can be analyzed by Eq. 14 to obtain a "quasiharmonic" entropy. Is this approach 
reliable? Some optimism was offered by Cowley that quasiharmonic theory would be 
useful for thermodynamic entropy of mildly anharmonic solids [324]. Wallace [46] 
and others have pointed out that the first order correction for the entropy from 
anharmonic frequency shifts is 

^ = -E4)^P. (166) 

kj 

where A(kj,T) is the energy shift of the phonon kj, and n £ (T) in the Planck occu- 
pancy factor. Using Eq. 14 for the phonon entropy, with the frequencies obtained 
from inelastic scattering measurements at elevated temperatures, is equivalent in 
first-order to using a harmonic entropy from the low-temperature phonon DOS plus 
the contribution from Eq. 166. 

This approach does not address any modification to the entropy by effects of 
phonon damping. Since the damping scales with the cubic term in the potential, 
by symmetry there is approximate cancellation of contributions in phase space 
of the position coordinates of opposite signs. Classical damping causes only small 
shifts of oscillator energies. Damping does cause a change in phonon linewidth, and 
the broadening of a phonon energy about a mean value is a concern for accurate 
determinations of vibrational entropy. If g(e) is broadened by a damped harmonic 
oscillator spectral intensity, for example, the vibrational entropy calculated with Eq. 
14 will be altered. Nevertheless, if the broadening functions for two materials are the 
same, e.g., they have the same value of Q in Eq. 139, their difference in vibrational 
entropy will be unchanged in the classical limit. Also, for small broadening at least, 
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the broadening is approximately symmetrical in e, and the thermal factors are 
approximately symmetrical. 

Although not fully justified, it is typical to assess the full vibrational entropy at 
elevated temperature by using Eq. 14 with the g(e) measured at the temperature of 
interest. In the few cases tested to date, this approach seems to work well. Eriksson, 
Wills and Wallace [334] used neutron scattering phonon measurements at elevated 
temperatures by Stassis, Petry, Heiming, and co-workers [337, 338, 339, 340] on bcc 
and hep Ti, and bcc and hep Zr to calculate the total phonon entropy at temperature 
(using Eq. 14). This was subtracted from the entropy obtained from calorimetric 
measurements, giving the electronic entropy, which was compared to calculation. 
The ultimate agreement was excellent. For all but bcc Ti, the discrepancy at any 
temperature was of order 0.01 fee/atom, which is probably better than can be ex- 
pected from the experimental measurements. For bcc Ti at elevated temperatures, 
the discrepancy was 0.1 k-Q /atom or less. Both hep Ti and hep Zr have modestly 
large anharmonic entropies at their highest temperatures (0.35 and 0.43 k-Q /atom), 
so the accuracy of the vibrational entropy, obtained from phonon measurements at 
temperature and analyzed with the harmonic formalism, is most encouraging. 

For all materials, phonons make the largest contribution to the heat capacity 
at high temperatures, and the phonon DOS must be obtained accurately. In cases 
where the phonon DOS has been measured at elevated temperatures, and the ther- 
modyanamic properties assessed, the anharmonic part of the phonon entropy was 
generally in good agreement with that deduced from thermodynamic data. The cases 
of elemental Ni and Al [341, 342] were in good agreement with the tabulations of 
Wallace [50] , although the contribution from the magnetic entropy of Ni complicates 
the problem somewhat [343]. The agreement was also good for bcc V [48, 344], but 
what had been assigned to anharmonic entropy from thermodynamic measurements 
proved to originate from electron-phonon coupling. It is not easy to separate the 
effects of anharmonicity from the effects of electron-phonon interactions by exper- 
iment alone, but conventional wisdom has been that electron-phonon coupling is 
negligible at elevated temperatures. 

19.5. Anharmonic Trends of fcc and bcc Metals 

The temperature dependences of phonon frequency shifts and broadenings give in- 
sight into the cause of anharmonicity. The quality factor Q of a damped classical 
oscillator has no dependence on the amplitude of oscillation, so there is no tem- 
perature dependence to its lineshape. This proves inaccurate in general - phonon 
frequencies usually broaden with temperature. A damping force that depends on 
the velocity of the nuclei is not fundamental to phonon anharmonicity, of course. 
From quantum mechanical perturbation theory we expect the first-order lineshape 
broadening to originate with the cubic term in the potential, and this is calculated 
to be linear in T. Such behavior is sometimes reported, but at higher temperatures 
the broadening often increases faster than T 1 (e.g., [333]). Perhaps this can be 
understood as follows. The potential energy in a nearly-harmonic oscillator increases 
at l/2kx 2 ~ k^T, so the oscillator displacement increases as x ~ VT. The limit to 
this displacement is approximately x/a = 0.15, where a is the interatomic distance, 
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Interatomic Separation (rescaled) 

Figure 39. "Universal" interatomic potential of Rose, et al. [345, 346] (solid curve). Points were 
obtained from Mo, K, Sm, Ba, Cu. After [346]. 



because larger displacements lead to melting at T m . The maximum value of the 

cubic term in the potential energy therefore increases as 

It was proposed in the early 1980's that the interatomic potential of different 
metals, scaled by binding energy, bulk modulus, and with minima centered about 
the interatomic distance, would follow a universal curve [345, 346]. The model 
seemed successful, as shown in Fig. 39. If true, the ratio of the cubic term to the 
quadratic term in the interatomic potential should be the same for metals. 51 The 
constants in Eq. 158 that depend on the interatomic potential are k and b. Assuming 
that the effective b for coupled anharmonic quantum oscillators causes a constant 
energy spread at a given temperature, then from from Eq. 159 we might expect that 
Q/ y/m ~ ^ , or 

Q / T \ 3//2 

* oc(=-) (167) 



for all metals if the phonon linewidth arises from the effects of the cubic term in a 
universal interatomic potential. 

The phonon DOSs of several fee metals were measured at elevated temperatures 
by Kresch, et al. [341, 342, 343]. For fee metals, the high temperature DOS curves 
could be fit well by modifying the low-temperature phonon DOS curves in two ways. 
The thermal shifts of phonon energies were small, and modeled successfully as a 
simple rescaling of the energy axis. The phonon energy broadening was accounted 
for by a damped harmonic oscillator function with the same value of Q for all 
phonons. These simple approximations were surprisingly successful for fee metals 
such as Al, as shown in Fig. 40 [343]. The energy rescaling was small and not 



51 When this potential of Rose, et. al., is fit to a polynomial form, the coefficients are comparable 
for the quadratic and cubic terms in displacement from equilibrium. The quartic term depends on 
the range of the fit, but is an order-of-magnitude smaller. 
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Figure 40. Phonon DOS of aluminum at temperatures as labeled. The markers show the exper- 
imental DOS. The lines are best fits to the elevated temperature DOS curves from rescaling the 
energy axis and convolution with the damped oscillator function Eq. 139. After [342]. 
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Figure 4.1, Experimental results from inelastic neutron scattering experiments at elevated temper- 
atures, from which Q was obtained by fitting Eq. 139 to experimental data from fee Al, Cu, Ni as 
shown in Fig. 40. T m is the melting temperature. After [343]. 



reliable from experimental measurement, but the trend in the phonon linewidth 
is presented in Fig. 41. There is a surprising agreement between the experimental 
phonon lifetime broadening for fee metals and the predictions of Eq. 167 for a simple 
model of a "universal" interatomic potential. The reader is cautioned that the data 
were obtained from only three fee metals, so there are not enough results on phonon 
measurements to address the vibrational entropy in a systematic way. Nevertheless, 
Fig. 41 gives optimism that the phonon anharmonicity in fee metals may be possible 
to understand in an general way. 

For bec metals, the trends are less simple because the different phonons do not 
shift with temperature in the same way [343]. This is particularly true of the low 
transverse acoustic phonons, which are intrinsically soft for the bee structure, but 
vary considerably for different metals. The low transverse acoustic modes tend to 
soften much more strongly with temperature than the longitudinal modes, and are 
associated with a number of structural phase transformations such as martensite. 
The 2NN distance in the bec structure is only about 15% longer than the INN 
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distance, so interatomic forces between 2NN pairs are expected to be important 
in the bec structure. This is generally true, based on tabulations of experimental 
interatomic force constants, for example [347]. These 2NN force constants can have 
a strong dependence on temperature, as appears to be the case for bec Cr [47], where 
the longitudinal force constants decrease by a factor of two over the temperature 
range from 300 to 1773 K, and the low transverse phonon branch shows a prominent 
softening over this temperature range. 

19.6. Thermoelectric Materials 

Thermoelectric materials need to transport electrons efficiently, but conduct heat 
poorly. Suppressing thermal conductivity is an important design strategy for new 
materials, and anharmonic effects are believed to be a viable approach to do so. 
Skutterudite antimonides of typical composition Lai(Fe,Co)4Sbi2 have a cage of 
Sb atoms around sites filled by Fe,Co, and around other sites filled by La. The 
addition of the caged La atoms suppresses substantially the thermal conductivity, 
and makes this family of materials interesting for service as a thermoelectric ma- 
terial. The localized vibrational modes of the La atoms were found by inelastic 
neutron scattering (using neutron- weighting to advantage in the experiment), and 
their vibrational levels were correlated to the heat capacity of the material [348]. 
The existence of resonance modes of the La atoms (sometimes called "rattling" ) is 
a characteristic of these materials, but with harmonic vibrations it is not obvious 
how thermal conductivity is suppressed. Subsequent work indicated that although 
the motions of the La atoms are mainly harmonic, with larger vibrational amplitude 
they couple more strongly to the acoustic modes responsible for thermal transport 
[349]. A possibility is that the local modes of the rattling atoms hybridize with the 
acoustic modes at small wavevector, k, effectively reducing the range in k of the 
linear dispersion of acoustic modes. This would suppress the number of acoustic 
modes available to transport heat. Another thermoelectric material, Zn4Sb3 seems 
to show local modes of dumbell rotation that suppress thermal conductivity [350]. 
These modes appear to be strongly anharmonic, and may scatter phonons out of 
the acoustic modes that transport heat. 



20. Electron-Phonon Interactions in Metals 

20.1. Fermi Surfaces 

20.1.1. Fermi Surfaces and Electron Screening 

Electron states near the Fermi level affect the interatomic force constants of metals. 
The electrons and holes screen charge disturbances, minimizing the electrostatic 
energy. The high effectiveness of this screening helps distinguish metals from insula- 
tors. Without screening by conduction electrons, the transfer of an electron between 
neighboring atoms of differing electronegativity would cause a Coulombic attractive 
force that falls with distance by about 1 eV/ A per A. This contributes about 16N/m 
to the interatomic force constant, which is comparable to the interatomic force 
constants themselves. Solutes with substantial differences in electronegativity would 
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then cause enormous distortions of interatomic distances and enormous changes in 
the phonon frequencies. The actual effects of alloying on lattice parameters and 
phonon DOS curves in metals are much smaller because screening is so effective. 
Since the atom displacements in a phonon distort the charge density, any changes 
in screening caused by alloying cause changes in the energy of a phonon. Likewise, 
it is important to look for changes in screening when considering how the phonons 
are affected by pressure or temperature. This Section 20 addresses how composition, 
temperature, and pressure alter the screening. 

Electron screening depends on the size and shape of the Fermi surface in fc-space. 
For example, the free electron model gives a spherical shape and an electron DOS 
that tends to be small at the Fermi level. Fermi surfaces are much more complex 
for transition metals, lanthanides, actinides, and their alloys [351]. Given the same 
topology of the Fermi surface, however, increasing the electron DOS at the Fermi 
level will increase screening, softening the phonons. Sometimes the simple parameter 
of the electron DOS at the Fermi level can be useful for understanding common 
features of how the phonon DOS is altered by effects of composition. The electron 
DOS at the Fermi level can also be used to understand the effects of temperature 
and pressure, although the effects of temperature require additional concepts. 

20.1.2. Fermi Surfaces and the Electron- Phonon Interaction 

To isolate the electron-phonon interaction, we reference the electronic energy to 
that of an undistorted, periodic structure without phonons and without thermal 
excitations of electrons. As discussed in Section 18.3, the "adiabatic" electron- 
phonon interaction requires simultaneous thermal excitations of both electrons and 
phonons. For a free electron gas with a spherical Fermi surface, we expect a strong 
contribution from of Eq. 154 when ki = 2/cf- Momentum conservation requires 
the three operators in Eq. 154 to give an annihilation of a phonon with wavevector 
of 2/cfz, for example, an electron moved into a state with wavevector —kpz, and 
an electron removed from a state with wavevector kpz. Both the electron states are 
at the Fermi level, so at finite temperature there is a good chance that the initial 
states are occupied and unoccupied as required. For a free electron gas, if ki differs 
significantly from 2kp, one or more of the electron states is very far from the Fermi 
surface, and excitation is improbable, even at elevated temperature. 

The Fermi-Dirac distribution blurs the Fermi surface. The thermal smearing 
function that accounts for thermal excitations across the Fermi surface is approxi- 
mately a Gaussian of standard deviation 2k^T. The available thermal energy also 
scales with T, so the energy of thermal electronic excitations scales as T 2 . The 
temperature-derivative of this energy is of course the electronic heat capacity, which 
is linear in T. These changes in electronic occupancies can alter the screening char- 
acteristics of the conduction electrons in some cases, but the effect can be small 
in others. For the spherical Fermi surface of a free electron gas, with increased 
temperature there is only a small change in the numbers of electrons and holes 
that participate in screening. It is often found that the elastic constants of metals 
have little explicit temperature dependence at a fixed volume. For many metals, 
most of the temperature-dependence of the elastic constants originates from the 
anharmonicity of phonon-phonon interactions, as described in Section 18.2, so effects 
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of "adiabatic electron-phonon interactions," as described in Section 20.3.1, are often 
neglected at high temperatures. 

Suppose the electron-phonon interaction is strong. At higher temperatures when 
more phonons are present, there is a proportionate increase in electron scattering 
and a shortening of the electron lifetime, r. A consequence is the lifetime broad- 
ening of the electron DOS by the characteristic energy AE = h/r. This lifetime 
broadening proves to be tricky to calculate because the increased scattering reduces 
the strength of the electron-phonon interaction that causes scattering. Nevertheless, 
in an approximate treatment [352] the resulting electronic broadening function is a 
Lorentzian function of full-width-at-half-maximum 2n\kBT. When A (the electron- 
phonon coupling parameter of Eq. 157) is approximately 1, this lifetime broadening 
exceeds the energy broadening caused by the Fermi-Dirac distribution. Furthermore, 
a Lorentzian function has slowly-decaying tails that sample a wide energy spread of 
the electron DOS. 

20.1.3. Kohn Anomalies and Electronic Topological Transitions 
For more complex Fermi surfaces, the situation is especially interesting when two 
sheets are parallel, or displaced by a common vector. This structure, sometimes 
called "nesting" of the sheets of the Fermi surface, allows many electron and hole 
states near the Fermi surface to be separated by the same wavevector. For a phonon 
of this wavevector, the screening of the charge disturbance is especially effective, 
reducing the energy cost of the phonon. This causes a reduction of the phonon 
energies in a narrow range of wavevector, called a "Kohn anomaly." Many Kohn 
anomalies have been reported (e.g., [353]). 

An "electronic topological transition" can occur when alloying fills (or depletes) 
the bands with electrons, and a feature of the Fermi surface is created or eliminated 
[354]. A feature of the Fermi surface for early transition metals is a hole pocket at the 
center of the Brillouin zone (the 'T-point," located at k = 0). Filling it causes a loss 
of states at the Fermi level and a loss of electrons that can participate in screening. 52 
There are a number of discussions about how electronic topological transitions alter 
physical properties of metals and alloys such as elastic constants [355, 356, 357, 358, 
359, 360]. Electronic topological transitions alter phonons by changing the adiabatic 
EPI. There have been relatively few systematic experimental studies of these effects, 
and most previous work was focused on selected phonon branches [361, 362, 363, 
364]. An early study of phonon dispersions in single crystals of bec Nb-Mo alloys 
showed large changes in phonons at a composition of approximately 66 at.% Mo 
that were interpreted with a rigid band model as an electronic topological transition 
caused by filling the hole pocket at the T-point [361]. At higher Nb concentrations 
the picture is more complicated, since an electron "jack" appears in the Fermi surface 
of Mo around the T-point. Adding Nb to pure Mo causes some phonons to stiffen 
and some to soften [361]. More recently an electronic topological transition at 30 
at.% Mo was calculated from first principles methods, but in the virtual crystal 
approximation, and was argued to account for some of the phonon anomalies near 



52 The electron DOS of vanadium, discussed below in the context of Fig. 43a, shows a drop of a 
factor of two at 0.35 eV above the Fermi level caused by filling the hole pocket around the T point. 
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this composition [365]. A similar effect is expected for W-Re alloys, but effects on 
phonons have not been studied [366]. 

Pressure can drive electronic topological transitions, and pressure is particularly 
convenient for computational work. A computational study on bcc Nb under pressure 
showed that topological transitions occur at 5 and 60 GPa [367], consistent with 
measurements of pressure-induced changes in the superconducting transition tem- 
perature [368]. A computational study on hep Zn showed an electronic topological 
transition at high compression that was responsible for a change in compressibility 
and the axial ratio of the structure [369]. Experimental evidence for a change in 
phonon frequencies in Zn was found by inelastic neutron scattering at a pressure of 
approximately 6.8 GPa [362], but the effect on elastic constants has been questioned 
[370]. An enhanced EPI after an electronic topological transition under pressure was 
recently reported in theoretical work on YB$ [371]. 

20.2. Phonons at High-Temperatures in Zirconium, Hafnium, Niobium, 
Molybdenum, and Chromium 

In a German-French- Austrian collaboration spanning several years, phonon dis- 
persions were measured for single crystals of elemental bcc transition metals at 
elevated temperatures [47, 339, 340, 372, 373]. Some of these results were discussed 
in Section 19.3 as they pertain to phonon-phonon interactions. Other aspects of this 
work focused on issues of lattice dynamics relevant to self-diffusion and soft-mode 
transitions. Some important thermodynamic results were also reported in this series 
of papers. In particular, the phonon DOS curves of Nb and Zr showed anomalous 
behavior with temperature - a slight stiffening of the phonons in both cases, most 
evident at the Van Hove singularity of the low transverse branch. Electron-phonon 
interactions were suggested as the source of this behavior. The metal Cr showed 
a very rapid softening of the phonon DOS with temperature [47]. An analysis of 
the lattice dynamics indicated that this originates from the thermal softening of the 
2NN interatomic forces, which are quite large for Cr at low temperatures. These 2NN 
forces have strong effects on the transverse modes of the bcc structure, and seem 
to be involved in the mechanism of shear transformations. The electronic origin of 
this softening in 2NN force constants is less clear, but it may be related to magnetic 
disorder [374]. 

A computational study of phonon-phonon interactions in bcc Mo [104] showed 
some discrepancies from the temperature dependence of phonons measured at the 
H point of the Brillouin zone [102]. A slight stiffening with temperature of phonons 
at the H point was observed experimentally, but the computational work (without 
electron-phonon interactions) predicted a large softening. The authors attributed 
this discrepancy to effects of Fermi surface nesting. 

20.3. Phonons in Vanadium and Its Alloys 

20.3.1. Temperature Anomalies of the Phonon DOS in Vanadium 
The phonons in pure bcc vanadium have an anomalous temperature dependence - 
Fig. 42 shows that there is no phonon softening to temperatures as high as 1000 
K [48], contrary to predictions of quasiharmonic theory (Section 6.3). On the other 
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hand, the phonon DOSs of V-6at.%Pt and V-Co show a normal softening with 
temperature that is fully consistent with the quasiharmonic entropy of Eq. 83. It is 
interesting that Nb solutes have a minimal effect on the anomalous lack of phonon 
softening with temperature, but the elements Nb and V are isoelectronic, and bcc 
Nb and bcc V have very similar Fermi surfaces. 

Figure 43a shows the electron DOS of V and three alloys, calculated by first- 
principles DFT methods [344, 355]. The addition of Nb has a negligible effect on 
the electron DOS, whereas Pt and Co both suppress the electron DOS at the Fermi 
level. The reduction in the number of states at the Fermi level should suppress the 
screening of charge density disturbances that occurs in the thermal motions of atoms, 
thus causing a stiffening of the phonons. This can be seen in the low-temperature 
data of Fig. 42, but is discussed extensively below. 

Figure 43b shows the thermal effects on the electron DOS of V and V-Co [344]. 
The broadening from the Fermi-Dirac distribution is well known, as is the temperature- 
dependence of the electron chemical potential, but the broadening from the electron 
scattering by phonons is less easy to calculate, and its curve at the bottom of the 
figure must be considered approximate [352]. Nevertheless, both these effects are 
large, and when convoluted with the electron DOS there is a considerable smoothing 
of the electron DOS and substantial changes in electron density near the Fermi 
surface. In essence, the effect of temperature on the electron DOS of vanadium is to 
reduce the electron DOS at the Fermi level in a way that is similar to the effect of 
alloying with Co. Both result in fewer electrons that can participate in the screening, 
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Figure 43. (a) Electronic DOS for V and V15X1 computed from first principles with VASP. (b) 
Electronic DOS for V and V15C01 at K and with broadening effect induced by the EP1 at 1000 
K and shift in ^{T). Also shown are the derivatives of the Fermi-Dirac function ~df/dE and the 
Lorentzian representing the EPI broadening for V, both for 1000 K. After [344]. 

causing a stiffening of the interatomic force constants. 53 For the V-Co, Fig. 43b shows 
less change of the electronic DOS with temperature, and the phonon DOS of V-Co 
in Fig. 42 shows normal quasiharmonic softening. For pure V, however, the loss of 
screening with increasing temperature suppresses the softening, effectively stiffening 
the phonons. The effect is approximately as large as the effects of quasiharmonic 
softening, although of opposite sign. The two effects cancel approximately, giving a 
phonon DOS that is largely unchanged over a wide range of temperature. 

It was indeed surprising that the effects of electron-phonon coupling in vanadium 
extend to temperatures as high as 1000 K, since electron-phonon coupling is typically 
considered in light of phenomena such as superconductivity, and the electron-phonon 
interaction is typically argued to be negligible at elevated temperatures [50]. 

20.3.2. Large Phonon Stiffening in Vanadium Alloys 

Inelastic neutron scattering experiments were used to measure the phonon DOS 
for bcc solid solutions of V with 3d, 4d, and 5d transition metal solutes [377, 355]. 
Results, organized by columns in the periodic table of elements, are presented in 
Fig. 44. The middle curves, from solutes in the V column of the periodic table, show 

53 From prior analysis of the electron-phonon interaction in vanadium, motivated by understand- 
ing its superconductivity, it seems that there are only modest differences in screening efficiency of 
the different electrons at the Fermi surface [375, 376]. 
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little change with alloying. The lower curves, from solutes one position to the left of 
V on the periodic table, show a slight softening of phonons with alloying. The top 
curves, from the elements in the Ni column of the periodic table, show a substantial 
stiffening. A vibrational entropy of mixing, AS 1 ™^, can be defined for a V-X alloy 
as the difference in vibrational entropy of the alloy and the vibrational entropy of 
the pure elements V and X, weighted by their atomic fractions. For dissolving 6.25 
at.% Pt in V, AS 1 ™^ is —0.272 &b /atom from Eq. 14. It exceeds the configurational 
entropy of mixing S , mix (0.0625) = +0.234 A; B /atom of Eq. 11. For a V-6.25% Pt 
alloy, the total entropy of mixing is therefore negative. Assuming the effects of solute 
atoms on the phonon spectrum are linear in the solute concentration, as expected 
at low concentrations, the entropy of mixing will be negative for all concentrations 
exceeding about 3.4% Pt. The configurational entropy of mixing will dominate at 
lower concentrations, owing to the logarithmic singularity of the first term in Eq. 11. 
Nevertheless, the composition where the total entropy of mixing is zero decreases 
approximately exponentially with the vibrational entropy contribution of the solute 
atom. A hypothetical solute causing twice the vibrational effect of Pt would have a 
negative entropy of mixing for compositions as low as 0.04 at.%. 

To perform the neutron weight corrections needed for thermodynamic accuracy 
of these results, a Born - von Karman model of lattice dynamics was developed for a 
16-atom cubic supercell containing 8 BCC unit cells of V atoms, except the central 
atom was a solute atom. The supercell dynamics allowed calculation of both the true 
phonon DOS and the experimental neutron-weighted phonon DOS (Sect. 8.3) by 
iterating interatomic force constants ^aa'iKi'n' to find a best fit to the experimental 
phonon DOS. For the alloy V-6.25at.%Pt, where the phonons stiffened considerably 
with respect to vanadium, a large Pt-V INN longitudinal force constant was found 
responsible for the change in phonon DOS of Fig. 44. The changes in the total DOS 
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were dominated by the partial DOS of V atoms having solutes in their INN shell 
(Fig. 45). 54 

Figure 46 shows a good correlation between the vibrational entropy of alloying 
and the electronegativity of transition metal solutes across the 3d row, and with Ad 
and 5d solutes down several columns of the periodic table. First-principles calcula- 
tions on supercells matching the experimental compositions predicted a systematic 
charge redistribution in the nearest neighbor shell around the solute atoms. The 
charge transfers between atoms were evaluated with the Bader algorithm (which 
tends to overestimate the charge transfer in metals, although it follows the Pauling 
electronegativity scale). Although this correlation points to the importance of elec- 
tronic rearrangements on phonon frequencies (as opposed to atomic size or mass, 
for example) , more fundamentally it is the screening of ion motions by the electron 
density that governs the restoring forces on the ions and hence the phonon frequen- 
cies (Sect. 20.1.1 and [379]). The d states of the late transition metal solutes Ni, 
Pd, and Pt are lower in energy than the rest of the alloy valence bands, leading to 
a split band for Pd and Pt solutes and considerable distortions of the band shape. 
Correspondingly, there is a systematic depletion of the number of states at the Fermi 
level for solutes toward the end of the transition metal series (about 20% for 6-7% 
solutes from the Ni column). This will reduce the number of weakly localized states 
available for screening the displacements of vanadium ions, causing higher energies 
for these displacements and higher phonon frequencies. The efficacy of screening 
depends on the available electronic excitations. The change with temperature of the 
screening can be attributed to an adiabatic electron-phonon interaction (EPI) as 
discussed next, even though the temperatures are quite high. 55 

54 Note also the resonance modes (Sect. 3.1) in the partial DOS curves of Pd and Pt around 
12 meV. 

55 It has been argued that the EPI is negligible above about a quarter of the Debye temperature 
of vanadium [17]. 
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Figure 46. Neutron-weighted vibrational entropy of alloying, AS*| b , plotted versus difference in 
Pauling electronegativity Ax (solid symbols) between solute and V, and Bader charge transfer on 
solute AQx (open symbols). The solid (dashed) line is a linear fit of AS*\ b versus AQx- After [377]. 



20.4. The Adiabatic EPI in Superconductors 

Figure 47 shows changes with temperature of phonon DOS curves of three com- 
pounds with the A15 structure [380]. The compounds VsSi and VsGe are well-known 
superconductors (with strong electron-phonon coupling), but V3C0 is not. The 
compound V3C0 shows a normal phonon softening with temperature that is approx- 
imately consistent with predictions of quasiharmonic theory. On the other hand, 
the phonons in VsSi actually stiffen with increasing temperature! This stiffening 
is weaker for VsGe, which has a lower critical temperature for superconductivity. 
Figure 48 shows the temperature dependence of the average phonon frequency for 
VsSi and VsGe. The average phonon energy first increases with temperature, then 
decreases at high temperatures. The low-temperature behavior is attributable to 
anomalously soft phonons owing to a strong electron-phonon interaction (EPI). 
With increasing temperature the EPI is suppressed, and quasiharmonic behavior 
is expected at the highest temperatures. 

The origin of the thermal suppression of the electron-phonon interaction is in- 
dicated in Fig. 49b, where the electron DOS at the Fermi level is shown to be 
strongly reduced with temperature for VsSi. This is caused by a thermal broadening 
of the electron DOS from the Fermi-Dirac distribution function, the thermal shift 
of the electron chemical potential, and from the shortened electronic lifetime at 
elevated temperatures. Figure 49a shows that in contrast, V3C0 does not have much 
structure of its electron DOS at the Fermi level, and it is found that the changes 
with temperature are small. 

Figure 50 shows phonon DOS curves for VsSi calculated with ab-initio methods 
[380]. The middle curves of Fig. 50, calculated with the thermal broadening of elec- 
tronic levels, are generally consistent with the experimental observation of phonon 
stiffening, especially for the low-energy modes in which the motions of V atoms are 



116 



B. Pultz 





200 400 600 800 1000 1200 



Temperature (K) 

Figure 48. Average phonon energy of A15 compounds of V~3Si, V~3Ge, and V3C0 as functions 
of temperature (from inelastic neutron scattering measurements). Markers are results from INS 
measurements. Dash-dotted lines correspond to conventional quasiharmonic (QH) behavior. The 
QH curves were offset vertically to match the experimental data at the highest T measured. After 
[380]. 



dominant. In contrast, the top curves show that volume expansion alone causes a 
softening of the phonons corresponding approximately to a rescaling of the energy 
axis. This is as expected with the quasiharmonic model, but is inconsistent with the 
experimental trend. At the bottom of Fig. 50 are phonon DOS curves calculated 
with consideration of the effects of both thermal expansion and the adiabatic EPI. 
These agree very well with the experimental trends. It is interesting that the lifetime 
broadening function required to fit the experimental trend is a bit larger than ex- 
pected, which suggests that the electron-phonon coupling parameter A is suppressed 
with temperature. 
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Figure 49. (a) Electronic DOS for A15 V 3 Si, V 3 Ge, and V3C0 from first-principles VASP cal- 
culations, (b) Electronic DOS for V 3 Si broadened by EPI at 25 K, 300 K, and 1000 K, using 
electron-phonon coupling parameter A = 1.0 at all T. Energies are with respect to the electron 
chemical potential fi(T). Inset: broadening function at 300 K and 1000 K. 
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Figure 50. Top: Phonon DOS for VaSi calculated from first principles for two different volumes 
(Vhigh = 1-03 x Viow). Middle: Calculated DOS for different electronic broadenings a. Bottom: The 
same as above for two different (V and a) combinations (cri ow = 50meV and 0hi g h = 200 meV). 
The calculated DOS was obtained from the neutron-weighted contributions of V and Si species and 
convolved with the experimental instrument resolution. The markers are experimental results for 
V 3 Si at 25 and 1023 K. 



Perhaps the most surprising result is that the EPI is thermodynamically impor- 
tant to 1000 K in A15 VaSi and VsGe, defying much conventional wisdom and many 
calculations of the effect of the EPI on the free energy of metals [50, 352]. Such a 
possibility was suggested by a few early theoretical papers that addressed EPI effects 
on the heat capacity [381]. More recent calculations on free-electron metals without 
structure in the electron DOS near their Fermi levels found effects of the adiabatic 
EPI at high temperatures, but these effects were small [382, 383, 384] (discussed 
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in Section 20.5). Other studies of early transition metals [385] and superconductors 
[386] suggested the importance of electron-phonon interactions on phonon dynamics. 

Phonons in the C32 phase of MgB2 have been studied thoroughly at low temper- 
atures, since the (nonadiabatic) EPI is the origin of its conventional superconductiv- 
ity. It is possible to fully substitute Al for Mg [387], obtaining the same C32 structure 
across the (Mg,Al) 11 B2 "quasibinary" system, and C32 AIB2 is not superconducting. 
Phonons have been studied with different Al and Mg concentrations [388], and the 
superconducting Mg-rich alloys show considerable softening of the high-frequency 
optical modes, consistent with a strong EPI and larger vibrational entropy. Early 
work indicated that superconductivity in MgB2 [389] originated with a strong EPI 
caused by anharmonic E2 g phonons. Later this interpretation has been called into 
question, and it has been argued that MgB2 is not particularly anharmonic [390]. 
Furthermore, the adiabatic approximation itself has been questioned because the 
electrons may not respond instantly to nuclear displacements [390]. A limitation of 
the previous work on phonon anharmonicity in MgB2 is that the measurements were 
at temperatures of 300 K and below, such as the extensive Raman spectrometry work 
that has been performed on MgB 2 and related materials [388, 391, 392, 393, 394]. 

20.5. Effects of Pressure, Composition, and Temperature on the 
Adiabatic EPI and the Phonon DOS 

For understanding how changes in the adiabatic EPI cause changes in phonons, we 
need some level of information about changes in the electronic structure near the 
Fermi level. If similar changes are caused by either temperature, composition, or 
pressure, we expect similar changes in electron screening, and similar effects on the 
phonons. For a simple example, increasing the band filling by alloying should have 
the same effect on the elastic moduli as a band filling caused by a reduction in 
volume with pressure. 56 For simple free-electron metals the effect of temperature is 
different, however, being a smearing the electron DOS. 

Bock, Coffey, and Wallace considered the adiabatic and non-adiabatic electron- 
phonon interactions with results from a pseudopotential perturbation theory [382, 
383, 384]. In essence, free electron theory was extended with perturbations from 
pseudopotentials, and perturbations from nuclear displacements, giving effects from 
the EPI that were second order in the pseudopotential. The non-adiabatic EPI was 
found to be small, in the range of 10~ 5 to 10~ 4 eV (but increasingly larger for K, 
Na, Pb, and Al), and weakly decreasing with temperature. The adiabatic EPI was 
calculated with phonon dispersions from experimental lattice dynamics, and with 
simple free electron dispersions for the electronic structure. For a reference state 
of a perfectly periodic crystal with no electronic excitations, the dominant effect of 
temperature was from the Fermi-Dirac distribution, f(T), which gave a factor of 
T 2 , and from the phonon occupancy factor n(e,T), which gave a factor of T 1 at 
modest temperatures. The numerical results for the adiabatic EPI were negligible 
at low temperatures, but increased as T 3 and exceeded the non-adiabatic EPI at 
the melting temperatures. Except for Al, however, the non-adiabatic EPI was a 

56 For a free electron gas with a spherical and featureless Fermi surface, we need only the electron 
DOS at the Fermi level to predict, for example, that the bulk modulus scales with the electron DOS 
at the Fermi level. 
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small effect. For Al at melting, the adiabatic EPI was 1CT 3 eV, which corresponds 
to an entropic effect of order 10 -2 fee/atom. It seems that the adiabatic EPI for free 
electron metals is not so important thermodynamically, but it is within an order of 
magnitude of being so. 

Metals with sharp features at the Fermi surface have larger thermal effects from 
the adiabatic EPI. The thermal broadening of the electron DOS can be understood as 
taking sharp peaks near the Fermi level and leveling them, or alternatively as filling 
states at low points in the DOS. We might expect that the detailed 3D structure of 
the Fermi surface could be important for understanding how phonon frequencies are 
affected by the adiabatic EPI, but the thermal effects of Sects. 20.3.1, 20.3.2, and 
20.4 showed good correlations to the simpler density of states at the Fermi level. For 
example, for A15 VsSi Fig. 49b shows that the Fermi level is at a peak in the electron 
DOS, so the thermal smearing and the shift in electron chemical potential reduce 
markedly the numbers of electrons and holes available for screening, reducing the 
EPI. Semiquantitatively, at least, this contributes more to the stiffening of phonons 
than in the case of pure bec V, where Fig. 43 shows a more mild effect of temperature. 
The electronic DOS curves of the A15 compounds and bec V are similar, however, 
with a Fermi level crossing the upper edge of a peak in the DOS. 57 

Nevertheless, generalizations of the effects of pressure, composition and temper- 
ature may be risky. In particular, under pressure the electron states remain well 
defined and long-lived. With increased pressure a complex Fermi surface is expected 
to remain anisotropic, and retain its complex features. With temperature, however, 
the electrons scatter out of their ground states owing to interactions with phonons, 
giving a lifetime broadening of the electron energy levels. A similar effect is expected 
from the effects of composition (if solutes are treated as scattering electrons from 
the states of a pure element). Intricate structures of the Fermi surface at T = will 
be lost at high temperatures owing to thermal broadening of the electron energies, 
however. The lifetime broadening causes a more interesting effect on the states at 
the Fermi level that participate in the EPI: 

— The lifetime broadening of the states at the Fermi surface owing to adia- 
batic electron-phonon interactions is dependent on the strength of the electron- 
phonon interaction in different directions in /c-space. 

— Those directions in fc-space with the stronger electron-phonon interactions will 
broaden more with temperature, suppressing the electron-phonon interaction. 

— The suppression of the electron-phonon interaction along directions in fc-space 
where it is strong will tend to equalize the strength of the electron-phonon 
interaction in all directions. 

The thermal smearing of the electron states and the increased isotropy of the 
electron-phonon interaction will therefore wash out the detailed structure of the 
Fermi surface. For a Fermi surface with one sheet, its low-temperature structure 

57 It is quite possible for the opposite effect to occur. When the Fermi level lies in a valley between 
regions of high electron DOS, it is expected that temperature will increase the number of electrons 
available for screening, and cause a further softening beyond what is expected from quasiharmonic 
behavior. 
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will evolve towards a fuzzy sphere at high temperatures. For a Fermi surface with 
isolated pockets, the anisotropy may persist to higher temperatures, but the detailed 
systematics require further research. 

Owing to the complexity and different structures of Fermi surfaces, further studies 
are needed to understand how well the effects of temperature on vibrational entropy 
can be predicted with only the electron DOS near the Fermi level. Furthermore, 
some conceptual issues such as how best to account for the lifetime broadening 
of electron states remain a challenge. Nevertheless, our present understanding of 
the concepts suggests that studies of phonons at high temperatures may provide 
useful information about the EPI, and perhaps about superconductivity, although 
the latter is more sensitive to the anisotropies of the EPI. 
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21. Systematic Trends of Vibrational Entropy 

21.1. Multiple Effects 

For solid-solid phase transformations, vibrational entropy has several contributions 
that are responsible for systematic behaviors. 

— It is important to count the different types of atom pairs in an alloy, and 
how they change in the phase transformation. The interatomic forces for INN 
pairs must be known. In many cases the effects of the INN triangle and 2NN 
pairs are important, but there is cause for optimism that larger clusters make 
rapidly-diminishing contributions. 

— Bond stiffnesses depend on bond lengths, and these effects can be comparable 
to differences in stiffness between different atom pairs. Effects of bond angle are 
not fully explored, but may be important for materials with d- and /-electrons. 

— Comparisons of bond stiffnesses are more difficult when there is a change in the 
crystal structure. The 2NN force constants are much larger in the bec than the 
fee structure, for example. 

— With changes in composition, not all effects can be attributed to stiffness 
changes of local bonds, because changes in the electron density at the Fermi 
level cause non-local effects that depend on the overall composition. 

— Alloying and pressure can be understood as changing the electronic structure 
near the Fermi level, which causes changes in the screening of charge distur- 
bances that accompany the vibrations of atoms. Modern DFT calculations can 
predict these changes quantitatively in many cases. 

— "Quasiharmonic behavior," is a classic thermal response that can be predicted 
from the bulk modulus and coefficient of thermal expansion. It reflects how the 
bond stiffness depends on bond length, and originates with the non-quadratic 
parts of the interatomic potential. Nevertheless, the interatomic potential is 
approximated as a harmonic one that changes with temperature. 

— Anharmonicities from phonon-phonon interactions originate from an interatomic 
potential that is unchanged with temperature. The phonon lifetime broadening 
is dominated by the cubic part of the potential, and may have a common 
behavior for fee metals. The phonon frequency shift depends on both the cubic 
and quartic terms of the potential, and trends are less obvious to date. 

— The adiabatic EPI accounts for how thermal excitations of electrons in metals 
interact with thermal excitation of phonons. For free electrons, the strength 
of the adiabatic EPI scale as T 3 , but larger thermal effects with different 
temperature dependences occur when there is big structure of the electron DOS 
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near the Fermi level. Pressure and composition can cause similar effects on the 
electron screening. 

21.2. HUME-ROTHERY RULES AND VIBRATIONAL ENTROPY 

The long previous list of phenomena affecting vibrational entropy may discour- 
age simple rules of thumb, but simple rules are used for other quantities in alloy 
thermodynamics. Assessing differences in vibrational entropy by rules analogous 
to Hume-Rothery rules [153, 395, 396] for metallic alloys risks predictions that 
are no better than the Hume-Rotherey rules for enthalpy, and likely worse. The 
Hume-Rothery rules of alloy phase formation are based on atomic properties: 1) 
electronegativity, 2) metallic radius, and 3) electron-to-atom ratio, with the first 
two factors being the most important. For phonon entropy, we should add 4) atom 
mass, to this list. 

— Electronegativity. This is a measure of the tendency for electron transfer to an 
atom during chemical bonding, which provides information on the strength of 
the bond and its stiffness. In metals this is only a part of the picture, because 
conduction electrons are highly effective in screening the charge disturbances 
that occur with electron transfer, but changes in the screening are also expected 
to scale with differences in electronegativity. Screening depends on the Fermi 
surface, but we expect systematic trends with electronegativity if, for example, 
the band structure is little changed with alloying. 

— Metallic Radius. Inserting a larger solute atom will cause a local compression 
in the host crystal. Since interatomic forces generally become stiffer under 
compression, alloying with large atoms tends to reduce the vibrational entropy. 

— Electron-to-Atom Ratio. This accounts for non-local Fermi surface effects in 
metals. The e/a ratio is expected to be a useful parameter if the rigid band 
model is reasonably reliable, and could be used to predict changes in screening, 
for example. 

— Mass. There is a general relationship for an oscillator of lo = ^Jk/M, where 
the frequency, cj, scales as the inverse root of the mass, M, assuming that the 
spring constant, k, is in fact a constant. Across the periodic table the atomic 
mass varies from 1 to 238 and beyond, and for a matrix having atoms of mass 
Mo, and solutes having mass Ms, with respect to the host the phonon entropy 
of the alloy scales with ln(Mo/Ms) [206]. No effects of mass are expected when 
a phase transformation rearranges the same atoms into a different structure, 
however. 
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22. Future Directions 

22.1. Opportunities for Experiment 

Measurements of the phonon DOS have long been possible with Fermi chopper 
spectrometers, but the experimental situation has changed markedly with the advent 
of new and efficient facilities for inelastic neutron scattering. There is now enough 
intensity so that dispersive information from coherent scattering can be measured 
reliably over wide ranges of Q and E. For powder samples this Q-information is 
rich in detail, but it is averaged over different crystallographic directions. With 
computational support, the S(Q, E) from polycrystalline samples can progress to- 
wards the more detailed information obtained from single crystals. Nevertheless, 
single crystal studies are possible with chopper spectrometers. With measurements 
at multiple crystal orientations, these studies should prove an efficient way to obtain 
phonon dispersions. Improved rates of data acquisition now allow more systematic 
measurements of phonons at elevated temperatures and modest pressures of a few 
GPa. 

A triple axis spectrometer is highly accurate in its energy and momentum mea- 
surements, and when a full set of phonon dispersions are measured on single crystals 
of ordered compounds, the phonon DOS derived from a lattice dynamics model 
is highly reliable. This technique is not so reliable for disordered solid solutions, 
unfortunately, because the dispersions are inherently broadened, and this broaden- 
ing need not be symmetrical in energy. Nevertheless, the biggest challenge remains 
obtaining single crystals of the material of interest. On the horizon are high inten- 
sity beams of long wavelength neutrons, which should offer new opportunities for 
studying long wavelength and low energy excitations in mesoscopic systems. These 
cold neutrons may also provide better energy resolution for phonon studies using 
phonon annihilations. 

Phonon dispersions can now be measured on very small single crystals by inelastic 
x-ray scattering. In spite of flux limitations, these experiments are practical today 
for many chemical compositions. This method could make significant advances with 
modest investments in energy analyzer systems. Inelastic nuclear resonant scattering 
is an incoherent scattering technique for determining phonon partial DOS curves of 
isotopes such as 57 Fe and 119 Sn, and is reasonably mature. Both these synchrotron 
methods offer the opportunities for measuring phonon information under extreme 
conditions of pressure and temperature, and these opportunities are starting to be 
exploited. 

22.2. Opportunities for Theory 

At low temperatures it is possible to use electronic structure calculations with the 
local density approximation to calculate phonon frequencies with good accuracy, 
but a 1% accuracy remains a challenge for many systems. Nevertheless, present- 
day DFT calculations are invaluable for finding trends in vibrational entropy with 
composition and pressure, and will see more widespread use. Better methods for 
electronic structure calculations are under active development. Improvements in 
ab-initio calculations on magnetic systems are expected as the exchange energy is 
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handled by methods beyond the local density approximation, but the cooperative 
nature of the magnetic transition or the magnetic susceptibility remains a compu- 
tational challenge. Density functional theory is not on solid ground when electronic 
excitations are considered, but there are developments in time-dependent DFT 
methods, and GW calculations give a more rigorous correlation energy. Path integral 
quantum Monte Carlo methods are also under development. These methods should 
prove valuable when considering matter in extreme conditions [397, 398]. Many of 
these computational methods will benefit from the ongoing rapid improvements in 
computing capabilities. Quantum molecular dynamics offers opportunities to handle 
anharmonic forces in a fundamental way, and with better computing resources these 
ab-initio molecular dynamics studies will be practical for larger systems and longer 
times. 

At modest temperatures, the suppression of the electron DOS near the Fermi 
level owing to the effects of electron lifetime is not well understood. The adia- 
batic electron-phonon coupling can be obtained from band structure calculations 
on supercells with atoms in static positions characteristic of thermal displacements. 
Results from an ensemble of such calculations give a broadening of the electron 
band structure that increases with temperature. Better calculations on the adi- 
abatic electron-phonon interaction can be used for more detailed predictions of 
thermal effects along different crystallographic directions. It may be possible to 
do this self-consistently, so the directions with the stronger EPI are more sensitive 
to temperature. 

The dynamics of disordered, nonlinear, or anharmonic systems have a long his- 
tory. The entropy of anharmonic systems will likely remain a challenge for many 
years to come. In the interim, it seems possible that better molecular dynam- 
ics calculations, perhaps ab-initio calculations, will play an increasingly important 
role in sorting out the dynamics of anharmonic materials and their temperature 
dependences. 

22.3. Outlook 

Much of the subject of vibrational thermodynamics of materials is based on the 
statistical mechanics of harmonic oscillators. It is fortunate that the harmonic ap- 
proximation gives a good start in understanding basic trends, because there are so 
many useful and powerful mathematical results for harmonic vibrations. Deviations 
from harmonicity are important, however, and understanding them is essential. High 
temperature behavior requires more concepts, and more complexity of analysis, 
but high temperature behavior is especially relevant to phase transformations in 
materials science and engineering. 

Experimental facilities for inelastic scattering of neutrons and x-rays have made 
large advances in the past decade, and offer to make routine the acquisition of 
systematic data. Many of the barriers to high quality experimental results are now 
on the analysis side of the work, rather than on the instrumentation side, as has been 
the case historically. Calculations of vibrational entropy are fruitful today, and will 
get better as computational resources improve and new approaches become available. 
The strengths and weaknesses of experiment and theory can often be overcome by 
working with data and computation in the same computing environment, using 
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the strengths of theory to support weaknesses in experimental measurements, and 
vice- versa. 

Many important systematic trends on the vibrational thermodynamics of materi- 
als are now known. Nevertheless, relatively few materials have been studied to date, 
and there is much more to do. The situation has a parallel to the atomic structure 
of materials, which has been an active area of research for a hundred years. The 
vibrational motions of atoms in materials, and how they originate from structure 
and in turn influence the thermodynamic stability of the structure, will remain rich 
for discovery for years to come. 
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Debye-Griineisen approximation, 12 
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